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EVALUATION OF DEFINITE INTEGRALS 
BY SYMBOLIC MANIPULATION* 



Abstract 



A heuristic computer program for the evaluation of real 
definite integrals of elementary functions is described. This 
program, called WANDERER (WANg's DEfinite integRal EvaluatoR), 
evaluates many proper and improper integrals. The improper 
integrals may have a finite or infinite range of integration. 
Evaluation by contour integration and residue theory is among 
the methods used. A program called DELIMITER (DEfinitive LIMIT 
EvaluatoR) is used for the limit computations needed in evaluat- 
ing some definite integrals. DELIMITER is a heuristic program 
written for computing limits of real or complex analytic func- 
tions. For real functions of a real variable, one-sided as well 
as two-sided limits can be computed. WANDERER and DELIMITER have 
been implemented in the MACSYMA system, a symbolic and algebraic 
manipulation system being developed at Project MAC, MIT. A typi- 
cal problem in applied mathematics, namely asymptotic analysis of 
a definite integral, is solved using MACSYMA to demonstrate the 
usefulness of such a system and the facilities provided by 
WANDERER. 



*This report reproduces a thesis of the same title submitted to 
the Department of Mathematics, Massachusetts Institute of Tech- 
nology, in partial fulfillment of the requirements for the degree 
of Doctor of Philosophy, August 1971. 
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Integral s. 

The limit concept Is basic to mathematical analysis. 
Being able to compute limits automatically greatly Increases 
the potential of a symbol manipulation system In doing 
analytical mathematics. In fact the evaluation of definite 
Integrals Is heavily dependent on the limit process/ as Is 
the expansion of functtons In power series and many other 
mathematical problems. Although the computation of limits 
has been studied previously to some extent [11/ 13]/ we 
describe a limit program called DELIMITER (DEflnltlve LIMIT 
EvaluatoR) which Is more powerful than previous programs. 
It is discussed tn detail In chapter 2. 

The problem of computing Indefinite Integrals 
symbolically by computers has been Investigated rather 
thoroughly. First among the computer programs developed for 
this purpose was SAINT (Symbolic Automattc INTegrator) by 
Slagle In 1961. A more powerful program named SIN 
(Symbolic iNtegrator) [23] was developed In 1967 by Moses. 
Theoretical work In this area Include Richardson's 
undecldabl 1 I ty result for a certain classes of tntegrals 
[22] and RIsch's decision procedure for determining the 
existence of the Indefinite Integral of a member In the 
class of elementary functions [23/ 2k, 25/ 26]. There Is a 



10 

comprehensive review by Moses on the progress during the 
past ten years In this particular area of symbol manipula- 
tion [21], 

However, hardly any work has been done In the direction 
of definite Integration. Reported here Is a first definite 
Integration program called WANDERER (WANg's DEflnlte 
IntegRal EvaluatoR). Experiments with DELIMITER and 
WANDERER are the principal subjects of this thesis. 
WANDERER has been Implemented In MACSYMA and makes use of 
many facilities provided In It. Some of these facilities 
are: Input and two-dimensional output, simplification, 
solution of polynomial and systems of linear equations, 
canonical rational function simplification and the SIN [20] 
Integration program. 

Evaluation of definite Integrals can sometimes be as 
easy as computing the Indefinite Integral of the given 
Integrand then substituting In the limits of Integration. 
WANDERER computes Integrals of this kind by using those 
parts of SIN that have being Implemented In MACSYMA. 
However, many Interesting definite Integrals are not 
obtainable In this manner. In some cases they may be 
Improper Integrals or their Indefinite Integrals do not 
exist. In other cases It Is easier to evaluate the definite 
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Integral directly than to obtain the Indefinite Integral 
first. 

A few ltmlt and Integration problems are listed below 
to give an Indication as to the scope of WANDERER and 
DELIMITER. The results obtained by the programs can be 
found In chapters 2, h and 5. 

f X + A X + B J XdX 

J k 2 L 2 15/2 

"" X ♦ 10 X +9 (X + X + 1) 

'" SIN(X) dX r °° 

-S T 1/3 

2 1X2 / e X LOG(X) dX 

-°° e (X + 1> J° 

1 dX ^ dX 

2 L 2 1/2 

X - 3 3 X (X - 9) 

2 f 1/2 2 

COS (X) - SIN(X) dX / LOG (X) dX/ X 

'o Jl 

2 4 1/2 I* 1/X 

LIMIT X U X + 5) - 2X LIMIT (A X + 1) 
X-*+<» X-*0 

2 1/2 2 
X (X + 1) X 
LIMIT e - e 

X-+a> 
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Integrals with an Infinite range are discussed In 
chapter h, those with a finite range In chapter 5. 
Important routines and algorithms are collected In chapter 
6. 

An application of the MACSYMA system to an asymptotic 
analysis problem Is Illustrated In chapter 7. The purpose 
Is to show how such a symbol manipulation system can be used 
to solve non-trlvtal problems that may occur frequently In 
applied mathematics. This application demonstrates the 
usefulness of the many faclllttes provided by MACSYMA and 
WANDERER. 

A timing experiment has been conducted to check the 
performance of WANDERER and DELIMITER. The results are 
Included In appendix C. Because of the limited character 
set of most computer consoles, some special symbols are 
needed to denote the frequently used mathematical functions, 
constants and operators. For example, * and ** are used to 
denote multiplication and exponentiation respectively. 
Appendix D contains a list of notations we shall use. 

From Richardson's undecldabll t ty results [22], we have 
shown that the convergence of a class of Integrals of 
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elementary functions Is recursively undec Idable. The proof 
Is In appendix E. 

The remainder of this chapter serves as an Introduction 
to algebraic manipulation systems. If the reader Is 
familiar with such systems he may proceed directly to the 
next chapter. 

An algebraic manipulation system Is a collection of 
computer programs designed to facilitate the solution of 
mathematical problems. Such a system has the ability to 
handle both symbols and numbers. This capability to mani- 
pulate symbols mathematically Is what differentiates these 
systems from the various computer subroutines which 
specialize In numerical analysts. 

In ustng symbol manipulation systems such as these, one 
usually Interacts with It In a time-sharing environment via 
a typewriter-like console. Normally, data and commands are 
typed In by the user. Results from the computer are sent 
back to the console for display. This arrangement attempts 
to provide a user the ease and flexibility of the pencil and 
paper he Is so accustomed to, while permitting the computer 
to assist him In his algebra and formal deduction from one 
step to the next. If the manipulations Involved are non- 
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trivial, so much the better. Some of the manipulations such 
a system can provide are: GCD calculat Ions, factoring, ra- 
tional function arithmetic, matrix manipulation, solution of 
algebraic equations, solution of systems of linear equa- 
tions, Indefinite Integration of elementary functions. 

Let us take a closer look at MACSYMA as a 
representative model of other systems. MACSYMA receives 
Inputs In the form of linear character strings typed by the 
user. FORTRAN-llke notation Is used for the Input. For 
example, 

The user types: if*(X**2+X + l)/Y@ 

The @ sign signifies the end of a command string. As a 
result of this command, a two dimensional display Is 
returned. 

1* (X + X ♦ 1) 
MACSYMA types: 

Y 

Inside MACSYMA, expressions are represented by list 
structures In a prefix notation (common to many systems). 
For example the expression above would be represented 
Internally by something like 

(TIMES 1* (PLUS 1 X (EXPT X 2)) (EXPT Y -1)) 
Inside MACSYMA. 
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Other forms of Internal representations are possible. 
The choice of Internal representation Is an Important aspect 
of the design of a system. Of course more than one form of 
Internal representation can be employed In a single 
algebraic manipulation system. In fact MACSYMA has a 
special internal representation for rational functions which 
Is used to gain efficiency during certain polynomial and ra- 
tional function manipulations. 

One simple application of MACSYMA Is factoring the 
polynomial 

3 2 

P(X) = X + l» X -11 X - 30 

A user who wants to factor P(X) using MACSYMA types: 
FACT0R(X**3+U*X**2-ll*X-33)@ 
This input command causes the factorization of PCX) over the 
Integers. The output is 

(X + 2) (X - 3) (X + 5) 

While this problem may seem easy, factorization of 
polynomials of higher degree can be very difficult to do by 
hand. Indeed, algebraic manipulation systems can be most 
helpful when one wants to manipulate complicated functions 
and expressions. MACSYMA can carry out accurately with 
great speed: summation over Indices, expansion of products 
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and powers, calculation of large determinants, Inversion of 
large matrices. 

Applications In pure mathematics Include computations 
In: number theory, group theory, Lie algebras and set 
theory. 

Research In the field of symbol and algebraic manipula- 
tion has led to many new results. The fast GCD algorithm 
[k] and the finite field arithmetic polynomial factoring 
algorithm [2] are two examples. Continued work In this area 
will, hopefully, result In computer systems which are 
Increasingly valuable to engineers and mathematicians. 
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CHAPTER I ! 
DELIMITER 

0. Introduction 

The limit concept Is fundamental to mathematical analy- 
sis. Basic concepts such as the rules of differentiation 
are derived from limiting processes. More complex problems 
such as Improper Integration/ convergence of series/ series 
expansion of functions and contour Integral s, to name a few/ 
also require the computation of limits In their solution 
process. 

Therefore, one can expect that automatic computation of 
limits would greatly Increase the capability of a symbolic 
mathematics system In doing analytical mathematics. The 
programs described below provide such a capability In the 
MACSYMA system [16]. 

Automating the computation of limits has been studied 
previously to some extent. Fenlchel [11] discussed certain 
decidability problems of limits and provided/ In the FAMOUS 
system, some baste routines for computing two-sided limits. 
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In his thesis, Iturrtaga [13] worked on both one and two- 
sided limits. In addition to L'HospItal rule, he applied 
asymptotic analysis to quotients of polynomials essentially 
by replacing those polynomials by their leading terms. This 
can be done when the variable In the limit approaches Infi- 
nity so that the leading term Is asymptotically equivalent 
to the original polynomial. Moreover, he discussed limits 
of sequences. 

DELIMITER Is more powerful than either Fenlchel's or 
Iturrlaga's limit program. L'Hospltal's rule is a basic 
method used for Indeterminate forms. In addition, this pro- 
gram employs a fast routine for limits of rational func- 
tions. It has an efficient algorithm for a class of expres- 
sions called RP-expressIons which Involve radicals of poly- 
nomials. It also applies the method of reducing complicated 
expressions by replacing subexpressions with asymptotically 
equivalent expressions. The method of comparing orders of 
infinity of expressions and several other heurtstlc methods 
are used. In some cases, power series expansions are employ- 
ed to obtain the limits. 

Thts chapter Is based on a paper by the author [29] 
presented at The Second Symposium on Symbolic and Algebraic 
Manipulation, Los Angeles, March, 1971. 
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1. Definitions and Symbols 

DELIMITER Is designed for finding the limit of a 
single-valued function f(X) of a real or complex variable X 
as X approaches some limit point. One-sided limits can be 
computed If X and f(X) are real-valued. The classes of 
functions allowed by this program Include rational/ radical, 
logarithmic and exponential functions, and also the trigono- 
metric functions SIN, COS and TAN, and the hyperbolic 
functions SINH, COSH and TANH. The general form of a user 
command Is 

LI Ml T( express I on, va r I ab 1 e, value, direct! on )@ 

with the fourth argument optional. PLUS as a fourth argu- 
ment indicates the one-sided limit from above, MINUS from 
below. The absence of the fourth argument tndlcates no 
restriction on the dlrectton of approach. One-sided limits 
are not allowed for complex-valued functions or limit 
points. 

If f(X) Is not continuous at X=a, the two-sided limit, 
L I M I T C f (X),X,a), does not exist. In such a case the symbol 
UND Is returned by the program as an answer. 



20 
The following are some symbols used by this program. 

ZER0+ 0+ ZERO- 0- 

UND undefined 

IND Indefinite but bounded 

INFINITY the point at Infinity In the complex plane 
INF and MlNF wMl continue to Indicate respectively 
positive and negative Infinity. Here are some examples of 
how these symbols are used Internally. 

UMIT(C0$(X),X,%PI) = -1+ ZER0+ 

LIMIT(SIN(X),X,INF) = IND 

LIMIT(1/SIN(X),X,INF) = UND 

LIMITC1/CX-*I>,X,*I> = INFINITY 
The direction from which a limit point Is approached Is 
Important. By use of ZER0+ and ZERO-, simplification rules 
such as 

1/C-l + ZER0+ + 1) — ► INF 
and 

l/(-l + ZERO- + 1) — ► MlNF 
are possible In DELIMITER. 
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2. Basic Rules and Outline of Algorithm 

Some of the rules used In this program are the so-call- 
ed "trivial" ones for limits of continuous functions, 
namely, the limit of a sum Is the sum of the limits, etc. 
There Is a limit routine for each of the functions SIN, COS, 
TAN, SINN, COSH, TANH and LOG (base %E). For other func- 
tions the rule 

LIMlT(f(g(X)),X,L) = LIMIT(f(Y),Y,LIMIT(g(X),X,L)) 

Is used. 

For one-sided limits where the variable approaches a 
point other than 0, a change of variable Is made to bring 
the point to the origin. For example, the limit 

LIMIT(f(X),X,a,PLUS> 
Is converted to 

LlMlT(f(a+Y),Y,ZERO+) 

There Is, of course, a complete set of rules governing 
simplification of the new symbols. They are all of the fol- 
lowing nature 
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IND/INF = -1/INF = ZERO- 

- ZERO+ = ZERO- - INF = MINF 
1/ INFINITY = -1/ INFINITY = 

MINF**2 = INF (ZERO+)+ZERO- = 
LOG(C+ZERO+) = LOG(C) + ZERO+ 
We will now briefly outline the algorithm/ and then go 
Into the details of some of the component routines. Let us 
cons Ider 

LIMIT <E(X),X,L) 

where L Is any number or symbol Including INFINITY/ INF, 
ZER0+ and ZERO-. If L = MINF, It Is set to INF by a change 
of variable Y=-X. Upon receiving the arguments, namely 
E(X)/ X and L, DELIMITER checks whether L-INFINITY or E(X) 
Involves %\ . If either or both Is true then a global 
Indicator CPLX Is set to the value TRUE which Indicates that 
the given limit Is to be evaluated over the field of complex 
numbers. Otherwise CPLX Is set to the value FALSE which 
signifies a limit problem over the reals. If CPLX=TRUE/ all 
use of INF and MINF are replaced by INFINITY and the notion 
of approaching a limit point from one side Is no longer 
val Id. 



The following diagram Is a simplified overview of 
DELIMITER. 
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DELIMITER 
SIMPLE LIMIT- 



-»RETURN ANSWER 



RATIONAL FUNCTION >RETURN ANSWER 



RP- FORMS 



♦-RETURN ANSWER 




BASIC LIMIT L'HOSPITAL'S SERIES 
RULES RULE EXPANSION 



ORDER 
OF INFINITY 



Fig. 1 

As indicated In fig. 1/ DELIMITER has four stages. 
Results of simple limit problems are Immediately returned by 
the first stage. Simple limits Include two cases: 1) E Is 
Independent of X. 2) E(X) = X. In the second stage/ a fast 
algorithm Is used to obtain limits of rational functions. 
The third stage employs an algorithm for limits of RP-forms. 
If the first three stages do not produce an answer/ the 
given limit problem enters the fourth stage which contains a 
variety of methods Including the four principal ones shown 
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In the figure. In this stage the program will first try to 
extract the numerator N(X) and the denominator D(X) of E(X) 
so that E(X) - N(X)/D(X>. If both H(X) and D(X) are RP- 
forms the limit Is computed by an algorithm specially for 
these forms (see section 5 for definition and algorithm). 
Otherwise the rule 

11m E(X) » Urn N(X) /l Im D(X) 

Is applied. If both 1 Im N(X) and 1 Im D(X) are or Infinite 
we have an Indeterminate form. In this case L'Hospltal's 
rule will be applied to E(X) with one exception: when L Is 
INF and 1 Im N(X) and Urn D(X) are both Infinite and both 
N(X) and D(X) contain exponential functions of X which tend 
to INF as X approaches L, In this case the method of com- 
paring orders of Inflntty, to be described shortly, Is ap- 
plied. If DELIMITER can not find an answer Indication of 
failure will be returned. 

Sometimes the program needs to know the value range of 
a symbolic parameter In order to compute the limit. In such 
cases, the program will query the user at his console. 

It Is Important to note that these programs do not 
store any table of limits Thus every limit obtained Is a 
result of computation. 
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Other methods available include : change of variable, 
discontinuity tests and analysis of the behavior of a func- 
tion near a finite point. A flowchart in appendix A pre- 
sents the flow of control in a more complete manner. 
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If n = Q, the answer Is c; 

If n < 0, then If CLPX = TRUE, answer is 0, but if 

CPLX is FALSE then the answer is ZERO+ provided 

the sign of c*L**n Is + and ZERO- In case the sign 

of c*L**n is -; 

otherwise, (n > 0) the answer Is INFINITY If 

CPLX=TRUE, and the answer Is INF or MINF, 

depending on the sign of c*L**n, if CPLX=FALSE, 
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U. The Limit of a Quotient of RP-forms 

An RP-form Is defined as an expression obtainable by 

combining polynomials and positive rational number powers of 

polynomials using the operators ♦, -, *. Here are some 

examples of RP-forms, 

3 1/3 

f(X) =SQRTU X + 5 X) + (X + 1) +X + C 

2/3 
g(X> = A X SQRTCX+1) + X 

Note that the definition of RP-forms does not allow 

nested radicals of X. For Instance the following Is not an 

RP-form. 



2/3 
(X + SQRT(X)) 



Now consider the limit problem/ 

LIMIT(E(X),X,U, ECX) -N(X)/D(X) 

where N(X) and D(X)j i are RP-forms and at least one of them 
Is Irrational. Note that N(X) or D(X) can be a constant. 

Let us define the operators EXPO and COEF by 
EXP0(N(X)) 

=htghest exponent of X In N(X), 
C0EF(E,N(X)) 

=the coefficient of E In NCX). 
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For example, 

EXPO(g(X)> =7/6, 
COEF(X**(3/2) / f(X)) =2 

EXPO Is obtained by a simple special purpose routine; COEF 
Is computed by using the RATCOEF routine In the rational 
function package of MACSYMA [16]. 

Now the limit Is computed by the following algorithm. 

(1) If L -MINF, X Is replaced by -X, L by INF. 

(2) If L =INF, the following asymptotic analysis Is 
made. 

I) compute a =EXPO(N(X)), b=EXPO(D(X) ) . 

II) Let NKX) and Dl(X) be N(X) and D(X) with 
polynomials under radtcals replaced by their 
leading terms, respectively. Compute 

a'-EXPO(Nl), b'=EXP0(Dl). 

III) If either a^a 1 or b^b 1 , then It Is usually 
sufficient to replace each radical In E(X) by 
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D(X) before using L'Hospttal's Rule. 

If N(L)=0 and D(L)?«0, 

the answer Is If CPLX=TRUE / otherwise the be- 
havior of D(L)*N(X) at L Is examined to decide 
whether 0, ZER0+ or ZERO- Is the answer. 

If D(L)=0 and N(L)?*0, 

the answer Is INFINITY If CPLX =TRUE, otherwise 
the behavior of N(L)*D(X) at L Is examined to 
decide whether INF or MINF Is the answer. 
Otherwise/ the answer Is E(L). 

For large X the series expansion 

2 1 

SQRTCX + 1) = SQRTC— + 1) X 

2 
X 

1 1 
= X (1 + + . . .) 

2 4 
2 X 8 X 
1 1 
= X + + . . . 

2 X 3 
8 X 

is convergent. It Is clear that this expansion method can 
be applied to an arbitrary polynomial to a positive frac- 
tional power. This fact establishes the validity of step 
(2-liI). 
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5. Indeterminate Forms and L'Hospital's Rule 

L'HospItal's rule Is applied to an Indeterminate form 
E(X) of the form 0/0 or INFINITY/INFINITY. It first 
calculates N 1 (X) and D'(X), then, 

If D'(X)=Q, L'Hospital's rule falls, otherwise, simplify 
N'(X)/D'(X) and evaluate 

LIMIT (N'(X)/D'(X),X,L). 

The expression N'(X)/D'(X) can be more complex than 
E(X) and successive application of this rule may lead 
nowhere. Therefore the number of times this rule Is called 
successively Is counted and the relative complexity of 
N'(X)/D'(X) to E(X) Is tested to decide whether to continue 
this approach. Our criterion of complexity Is based on the 
number of distinct nonratlonal components of an expression. 
If this number grows for three consecutive times, the ap- 
plication of L'Hospital's rule Is halted. 

The Indeterminate form 0*1 NF I N I TY Is transformed to 
either 0/(1/ INFINITY) or I NF INITY/U/0), depending on which 
Is simpler, before applying L'Hospital's rule. Other In- 
determinate forms such as 1**INFINITY, INF**0 and 0**0 are 
handled by the logarithmic reduction : 
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F(X)**G(X) = %E**(G(X)*LOG(F(X))). 

Theoretically the Indeterminate form (INF - INF) can 
always be reduced to 0/0 by rewriting It as 

(1) (1/ INF - 1/INF)/(1/!NF*INF) / 

but this method often makes the expression much more 

complicated. It Is useful/ though/ for expressions 

Involving trigonometric functions as can be seen In 

1 1 

LIMITC /X/0,PLUS) = 0. 

X SIN(X) 

DELIMITER transforms the given expression In this example to 

SIN(X)-X 

X SIN(X)/ 
then applies L'Hospltal's Rule to obtain the answer 0. 
Therefore/ the method (I) Is used for expressions Involving 
trigonometric functions. Other types of expressions can be 
dealt with more readily by comparing degrees of Infinity of 
the subexpressions or by series expansion. 
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6. On the Sign of Infinity 

If a function f(X) Is discontinuous at X=a, another 
function g(f(X)) Is not necessarily discontinuous at a. 
For example 

1/X A 

LIMITUA X + 1) ,X,0)= %E 

Therefore/ discontinuity may be encountered In the course 
of computing a limit of a continuous function. Consider 
1/ECX) with E(a)=0 and the limit problem 

LlMlT(l/E(X),X,a), 

If CPLX = TRUE, the answer Is INFINITY. Otherwise, a rou- 
tine named BEHAVIOR Is used to analyze the behavior of E(X) 
near X=a. If E(X) approaches from above or below as X ap- 
proaches a, then the answer Is INF or MINF, respectively. 
Otherwise, the answer Is INFINITY. 

The BEHAVIOR routine uses differentiation to analyze 
the behavior of a function near a point by Investigating Its 
slope or concavity at the point. This routine Is also used 
by other programs. One of these Is the program which com- 
putes limits of the trigonometric function tangent. In case 
the argument of TAN approaches £PI/2, say, It Is Important 
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to know the direction of approach. 

7. Heuristics for Comparing Orders of Infinity 

Let f(X) and g(X) be two real-valued functions which 
become positively Infinite as X approaches IMF. We define 
the symbol >> by 

Definition : g(X) >> fCX) If 

UMIT(f(X)/g(X),X,INF> =0. 

Using this symbol we can make the following brief 

table. 

... %E**2E**X>>%E*X>>X**n>>log(X)>>logClog(X)) ... 

Many limit problems with the variable approaching INF 
can be solved very efficiently by using this concept. Some 
of these problems can be difficult to solve by other me- 
thods. For example, L'Mospltal rule falls to compute the 
1 imlt 

LIMIT((%E+1)**X**2/%E**X,X,INF>, 

while the answer Is obviously INF. We Incorporate this con- 
cept by a routine, STRENGTH. It can classify the order of 
Infinity of the argument according to the following rules. 
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STRENGTH(c)=(c) 
STRENGTH(LOG(R(X)))= LOG 
STRENGTH(A**B(X))» EXP 
STRENGTH(P(X)**(m/n))= m*d/n 

where c Is Independent of X, m and n are positive Integers, 
d the degree of the polynomial P(X). Note the STRENGTH of a 
constant Is denoted by that constant Inside parenthesis and 
the STRENGTH of a polynomial Is a constant. B(X) and R(X) 
can be any functions In X, except those which can cause the 
relevant argument of STRENGTH to be simplified Into one of 
the other three cases. We can assume A to be 2E for If A ts 
any other expression It Is always possible to change the 
base to &E. STRENGTH of a sum Is the maximum of the 
strengths of the terms In the sum. 

There Is a basic comparison routine which knows the re- 
lations between the strengths of functions, 

EXP>>a>>LOG>Xb) and c>>d If c>d 

where a, b, c, and d are constants. 

Let N(X) and D(X) be two products satisfying 

LIMIT(N(X),X,INF)=INF 
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and 

UMIT(D(X),X,INF) = INF 

and no factor In N(X) or D(X) Is a sum or a trigonometric 
function. To compare the relative order of Infinity of 
N(X) and D(X), the following method Is used. The answer 1 
Indicates N(X)>>D(X), -1 Indicates D(X)>>N(X), no 
decision. 

Algorithm COMPARE(N(X) ,D(X) ) : 

1) Remove from N(X) and D(X) any commop factors. 

2) Apply the STRENGTH routine to each factor of N<X) 
and obtain the maximum strength SN. 

3) Do the same to D(X) and obtain Its maximum strength 
SD. 

k) If SN>>SD the answer Is 1, If SD<<SN the answer is - 
1. 

Otherwl se / 

I) Let 

Nl(X)=product of all factors with strength SN In 

N(X). 

Dl(X)=product of all factors with strength SD In 
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D(X). 

i!) If SN=SD=some constant/ set 

SN=STRENGTH(N1) and SD=STRENGTH(D1) . 

If SN>>SD / the answer is 1, If SN<<SD the answer Is 
-1. If both NKX) and D1(X) are polynomials In X, 
the answer Is the leading coefficient of the 
polynomial NKX) - D1(X). Otherwise If the variable 
INDICATOR has the value 1, this means each of N(X) 
and D(X) Is an exponent of an exponent of an expo- 
nential expression and they come as a result of 
recursive calls to COMPARE. For example N(X) may 
come from %E**2E**N(X) . In this case, the algorithm 
goes as follows. First INDICATOR Is set to then 
E=UMIT(N(X)-D(X)/X/ INF) Is computed and the answer 
Is 1 If E>0, -1 If E<0, otherwise. 

If SN=SD=EXP, 

The exponents A(X) and B(X) of Ml(X) and D1(X) 
are computed. That Is N1(X) =?E**A(X) and D1(X) 
=?SE**B(X). If both A(X) and B(X) are polynomials In 
X/ the answer Is the leading coefficient of the 
polynomial A(X) - B(X) tf this polynomial Is not a 
constant/ If It ts. If A or B Is a sum/ A Is set 



39 

to LIMITCA(X)-B(X),X,INF) and the answer Is 1 If 
A=INF, -1 If A=MINF / otherwise. In case neither A 
nor B Is a sum/ If both A and B are of strength EXP, 
set the variable INDICATOR to 1. The answer ts ob- 
tained by evaluating C0MPARE(A1(X),B1(X)>, where 
A(X) = 2E**Al(x) and B(x) = %E**B1(X). 

If SN=SD=LOG, 

If Nl(X) Is LOG(fCX)) and Dl(X) Is L0G(g(X)) then 
set SH to STRENGTHCf(X))/ SD to STRENGTH(gCX) ) . Now 
the answer Is If both SM and SD are constants/ 1 
If SN>>SD/ -1 If SN«SD. The answer Is otherwise. 

We shall next discuss the Indeterminate form (INF-INF) 
and see how COMPARE can be used In such situations. Let 

tt 
F(X)= £ fi (X) 

i=i 
where 

L 1 M 1 T ( f . (X)/X/INF) -INF or MINF 

for 1=1/2/ . . .n. 

LIMIT(F(X)/X/ INF) becomes Indeterminate when there exists I 

and j such that 

LIMlT(fj (X)/X,INF) » INF 



and 



LlMlTCf- (X)/X/INF) = MINF 

J 



The method of orders of Infinity can sometimes be ap- 
plied to this problem. The algorithm used Is as follows. 

The program forms two lists LI and L2 of the f j ' s which 
go to INF and MINF respectively. Then from each of LI and 
L2 a member of maximum order of Infinity can be obtained. 
Let these be r(X) in LI and s(X) In L2. Now if r(X) >>s(X) 
the answer Is INF, If s(X)>>r(X) the answer Is MINF. Other- 
wise, the problem can be very difficult. However the pro- 
gram tries to evaluate 

INF*LIM!T(F(X)/s(X),X, INF), 

which sometimes produces an answer. 

Now we w I 1 ] follow the major steps of the solution of a 
problem In a more detailed manner. Let 



A(X) = %E 



X*SQRT(X + 1) 



B(X) » %E 
Consider the problem L I Ml T( AC X) -B( X) ,X, ! NF) . First the 
Indeterminate form (INF-INF) is encountered through the fol- 
lowing steps. 



IMIT(A(X)-B(X) / X / INF) 



LIMIT(A(X),X,INF)-LIMIT(B(X),X,INF) 

LIMIT(A(X),X,INF) =%E**INF = INF 

LIMIT(B(X),X,INF) *%E**INF = I NF 

INF-INF 
Then COMPARE Is applied to A(X) and B(X), only to 
arrive at no conclusion about the relative order of Infinity 
of A(X) and B(X). 

COMPARE(A(X),B(X)> 

STRENGTH(A(X))=EXP 

STRENGTH(B(X))=EXP 

2 2 

COMPARE(X*SQRT(X + 1),X ) 

2 

STRENGTH(X*SaRT(X + 1)) = 2 

2 
STRENGTHU ) = 2 

Now the problem Is converted to the following form and the 
answer Is INF. 

INF*(LIMIT(A(X)/B(X),X,INF)-1) 

2 2 
LIMIT(X*SQRT(X + 1) - X ,X,INF) = 1/2 
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INF*UE**(l/2)-l) 
INF 

8. Examples 

A number of ltmtt problems solved by DELIMITER are 
Included In this section. They are presented In the form of 
actual Inputs and outputs of the MACSYMA system. Lines 
labelled (Cl) are Input or command lines and (Dl) output 
1 I nes or answers. 

(Cl) LIMIT(X**L0G(1/X) / . X, I NF > Q 
(Dl) 

(C2) (COS(X)-l)/( %E**X**2 -1)(3 

COS(X) - 1 



(D2) 



2 
X 

%E - 1 



(C3) LIMIT (D2,X,0)(a 

(D3) 

(Cl*) (1+A*X)**(1/X)@ 

(D4) (A X + 1) 

(C5) LIMIT(Dt*,X,0)@ 

A 
(D5) %E 

(C6) X**2*( lj*X**U + 5 )**(!/ 2 )-2*X**4@ 



1 
2 
l./X 
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(C6) X**2*U*x**lt + 5)**(l/2)-2*X**l*@ 

2 k k 

(D6) X SQRTU X + 5) - 2 X 

(C7) LIMIT(D6,X,INF>@ 

5 
(D7) 

(C8) 1/X-1/SIN(X)@ 

1 1 
(D8) 

X SIN(X) 
(C9) LIMIT(D8,X,0,PLUS)@ 

(D9) 

(CIO) fiE**(X*CX**2 + l)**(l/2))-^E**X**2@ 

2 2 

X SQRT(X +1) X 

(DIO) %E - %E 

(Cll) LIMIT(D1,X,INF)@ 

(DID INF 

(C12) UE**X+X*LOG(X))/(U**(X**3 + l)**(l/2) + LOG(X**4+X + l>>@ 

X 
%E + X LOG(X) 
(D12) 

3 

k SQRT(X + 1) 

LOG(X + X + 1) + |E 

(C13) UMlT(Dl»,X,INF>g 

(D13) 

(Cl'O l/(X**3-6*X**2+ll*X-6)@ 



kk 



(D1U) 



3 2 

X -GX +11X-6 



(CIS) LIMIT(D1U,X,2,MINUS)@ 



(D15) 



NF 



(C16) (X*SQRT(X + 5) + l)/(SQRTU*X**3 + l) + X)@ 



(DIG) 



X SQRT(X + 5) + 1 



SQRTU X + 1) + X 



(C17) LIMIT(016,X,INF)@ 



(D17) 



1 

2 



(C18) TAN(X)/LOG(COS(X))@ 



(D18) 



TAHCX) 
.OG(COS(X)) 



(C19) LlMlT(D18 / X / ?PI/2 / MlNUS)@ 



(D19) 



Ml NF 



(C20) (Z-*l*>;;PI/2)*Z*(Z-2*2PI*?;l)/(SIN ! H(Z)-%l)i 



kS 



(D20) 



?< P I %\ 

2 %PI %\ ) (Z - -■ --} 

2 

SINH(Z) - ?! 



(C2I) L!M(DiFF(!)20 / n,Z^I*%PI/2)Q 



(D21) 



- 2 %P 
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CHAPTER I I I 
WANDERER-AN INTRODUCTION AND OUTLINE 

1. Introduction 

The present chapter and the three subsequent ones are 
devoted to a complete description of WANDERER/ a heuristic 
program for evaluating definite Integrals over a real range. 

WANDERER cannot solve all definite Integration 
problems/ but It can certainly obtain solutions to a large 
number of Integrals. Many examples presented are actual 
problems taken out of graduate text books [1/ 5], The basic 
approach, tn WANDERER/ for evaluating an integral Is by com- 
plex contour Integration and residue techniques. This Is a 
very general method particularly useful In evaluation of 
Infinite integrals. Some other methods available to 
WANDERER are: substitution/ Integration by parts, diffe- 
rentiation with respect to a parameter, pattern recognition 
and table look-up, f I ni te-to-t nf in I te conversion, Introduc- 
tion of a parameter and partition and transformation of the 
range of Integration. 
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The user's command In MACSYMA for definite Integration 
Is the key word DEFINT. DEFINT takes four arguments as In 

DEFINT(exp, var, a, b) 

where exp Is the Integrand, var the variable of Integration, 
a and b the lower and upper limits of Integration. There Is 
no restriction as to what types of Integrands are allowed as 
Input, as long as they remain elementary functions. 
However, WANDERER requires the Integrand to be finite In the 
Integration range except possibly at the end points a and b. 

A few special symbols are used In WANDERER for some 
branches of the logarithm function. They are listed In 
appendix D. 

2. Outline of Approach 

Before going Into the details of methods and algorithms 
for evaluation of the many types of definite Integrals which 
will be discussed In the next three chapters, a quick look 
at the whole picture with emphasis on the flow of control is 
In order. 

WANDERER Is a heuristic program which computes definite 
Integrals by trying to apply one or more of the methods or 
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algorithms built Into It. The clues as to which of the 
methods to apply are obtained by examining the range of 
Integration and the form of the integrand. 



Pig, 1 on the next page serves as a simple outline of 
the flow of control in WANDERER. 



S=success 
F=fal lure 



WANDERER 



simple 



ntegral s-»»return answer 
F * 



normal izatlon 
of Integrand and 
Integration range 
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methods for 
Infinite Integrals/ 
applied according to 
range of Integration 



return answer 



methods for 

finite Integrals/ 

Including proper and 

Improper integrals 



more general methods 

transformation/ ^return answer 

change of variable £ 

F 




divergence test 



Fig. 1 



As Indicated In fig. 1. WANDERER can be said to have 
five stages. In the first stage simple Integrals are com- 
puted which Include cases such as a constant Integrand/ an 
Integral with equal upper and lower limits of Integration 
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and / if the variable of integration Is X, Integral with 
Integrands which are polynomials In X or %E**X. The 
normalization stage will make sure / by transformation If 
necessary, that (1) the lower limit Is less than the upper 
limit, and (2) constant factors, if any. In the Integrand 
are removed to be multiplied Into the final answer of the 
integral . 

The thtrd stage comprises most of the methods that will 
be described in chapters h and 5 and Is the work-horse of 
the whole program. Many methods In this stage are grouped 
according to the integration range in which they are 
appropriate. Finite integrals of rational functions are 
transformed Into infinite Integrals which will be evaluated. 
For other finite Integrals, an effort Is first made to 
determine whether they are Improper Integrals. The absolute 
divergence of a finite improper integral is tested before 
any attempt at evaluation. If the third stage fails to 
obtain any results, the given problem enters the fourth 
stage which contains more general methods appropriate to 
various types of Integrands. A transformation or change of 
variable is often done In this stage. When WANDERER runs 
out of methods, the convergence of the given integral will 
be challenged. In case it is divergent, WANDERER will so 
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Indicate In the output. Otherwise, Indication of failure 
will be returned. It may be argued that tests for 
convergence should be conducted before evaluation, but this 
Is not necessary because WANDERER Is so designed that 
whenever it produces an answer, the given Integral Is 
convergent. That is to say each Individual method has Its 
own convergence and divergence conditions built In. 

The limit program described In the previous chapter Is 
used whenever a limit computation Is needed In WANDERER. It 
may often be required In changing the variable of Integra- 
tion, computing residues and testing for convergence or 
divergence. 

3. Outline of the Computation of Residues 

One of the most powerful methods tn evaluation of 
deftnite Integrals Is the use of contour Integration through 
residue calculations. The Importance of contour Integration 
and the residue theory In the sequel warrants a brief 
summary of relevant facts from complex analysis [5]. 

(i) The Cauchy Integral Theorem: Let D be a simply 
connected domain and let f(Z) be analytic tn D. Let C be a 
closed contour tn D. Then 
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I 



f(Z) dZ » 



A simply connected domain Is Intuitively an open set of 
points with no holes or cuts In Its Interior. For example 
the disk |Z| < Is simply connected while the annulus 1 < 
|Z| < 2 Is not. The Cauchy Integral theorem Is one of the 
most Important tools In complex analysis. It Is by use of 
this theorem that path of Integration can be deformed. More 
precisely/ If CI and C2 are two different curves In D lead- 
ing from the same starting point v to the same end point w, 
as shown In fig. 2, 




Fig. 2 
and f(Z) Is analytic In D, then 



f f(Z) dZ = f f 



f(Z) dZ. 



Therefore/ the contour CI can be deformed Into C2. The 
above Is an Immediate consequence of Cauchy's integral 
theorem which gives 
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L 



f(Z) dZ - / f(Z) dZ « 

(If) Deflnltton of Residue: Let f(Z) be analytic In the 
punctured neighborhood of a point p. Let Cp be a small 
positively oriented circle with center at p. Then the 
residue of f at p Is defined as 






f(Z) dZ 



The value of this Integral Is the coefficient of the 

term 

-1 
(Z - p) 

In the Laurent expansion of f about the point P/ so this 
coefficient can also be regarded as the residue of f at p. 

(Ill) The Restdue Theorem: Let C be a simple closed 
positively oriented contour. Let D be a simply connected 
domain containing C and Its Interior. Let Al/ A2/ . . . , An 
be points Inside C. Except for Isolated singularities at 
Al, A2/ . . ./ An/ let f(Z) be an analytic function In D. 
Then 
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r n 

/ f(Z) dZ = 1%\ £ (Residue of f at A j ) 

The algorithm for residue computation used In WANDERER 
will be described briefly here. Detailed discussion of this 
algorithm can be found In sect. 6-1. 

Let f(Z) be a function of a complex variable Z, 
analytic everywhere In a domain D except for a number of 
poles In D. Suppose F(Z) can be written In the form 

f(Z) » U(Z)/V(Z) 
such that U(Z) Is analytic In D. This means. that poles of 
f(Z) are zeroes of V(Z). Suppose p Is a pole of order m of 
f(Z)/ the residue of f(Z) at p Is computed by the following 
al gorl thm. 

RESIDUE ALGORITHM: 



If m - 1/ compute as answer 

U(p)/V'(p) 
otherwise/ If V Is a polynomial/ 

(I) Set V to the quotient of V(Z)/(Z - p)**m 
which Is computed by long division. 
(ID Return the result computed from 

m-1 



(m - 1)! dZj 



U(Z) 

"vTzT 



Z=p 
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otherwlse / (m > I, V not a polynomial) compute and 
return as answer the limit 

m-1 



Urn 



d \ m 

— (z - p) f(z> 



(m - 1)1 \dZ/ 



DELIMITER Is used In obtaining such a limit. Methods 
for Infinite Integrals are discussed In the next chapter. 
Those for finite Integrals are Included In chapter 5. 
Important algorithms and algorithms common to many methods 
are detailed In chapter 6. 
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CHAPTER IV 
SYMBOLIC COMPUTATION 
OF INFINITE INTEGRALS 

0. Introduction 

The aim of this chapter Is to give a detailed account 
of the methods employed by WANDERER for Improper Integrals 
with an Infinite range. Contour Integration and residue 
computation play a very Important role In many of the 
methods. The algorithm for Integrals with rational 
Integrands/ which will be discussed first/ Is the most com- 
plete. Sections are formed according to the function types 
of the Integrand and ordered roughly In Increasing com- 
plexity. Quite a few examples are Included/ some of them 
with references to books or Integral tables Indicated. 

1. Infinite Integral of a Rational Function 

In this section methods for evaluating Infinite Inte- 
grals of a rational function will be presented. The range 
of Integration Is from either or minus Infinity to 
Infinity. Other possible Infinite ranges such as (a INF)/ 
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(MINF a), where a Is finite, can be converted to (0 INF) by 
a change of variable In the Integral. 
First Vet us consider 

L =/ R(X) dX, R(X) = P(X)/Q(X) 
J-co 

where P(X) and Q(X) are polynomials In X over the field of 
complex numbers. 

WANDERER requires deg(Q(X> )-deg(P(X) ) S 2 to Insure the 
convergence of the given Integral. If R(Z) has no real 
poles, then The Integral L can be computed by evaluating the 
contour Integral 




around a familiar seml-c I rcul ar contour In the upper complex 
Z-plane (fig. 1). One can easily prove that L = J as r 
tends to Infinity. 




*-x 



UfA. 
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Therefore 

(1) L = 2 ?;PI %\( £ Res(R) over poles of R Inside C) 

The residue subroutine outlined In the previous chapter Is 
used to compute the residues needed In this formula. If 
R(Z) has poles on the real axis then L Is divergent. 
However, If the real poles of R are all simple then the 
Cauchy Principal Value of the Integral L exists, and can be 
obtained by Indenting the contour C at these singular points 
on the real axis. For such Integrands WANDERER will compute 
this principal value, which Is given by 

(2) (P)L = 2 %P\ %\( £ Res(R) over poles of R Inside C) 

+ %P\ %\( Res(R) over simple real poles of R), 

as an answer to the Integral L. (P)L stands for the 
principal value of L. Whenever the answer Is a principal 
value, the message PRINCIPAL will be sent to the user first. 

Here is an example solved by WANDERER 
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(CI) (X**2+A*X + B)/(X**l* + 10*X**2 + 9)@ 

2 
X + A X + B 
(Dl) 

l» 2 
X + 10 X +9 

(C2) DEFINT(D1,X,MINF,INF)@ 

*PI B + 3 %PI 



(D2) ■ 

12 



The expression Dl has four simple poles: X = %\ , -%\ , 
3 %\ and -3 %\ . Using formula (1), only residues at X=%l 
and X=3 %\ need be computed. The residues are 

at X=%l (B-l+A *l)/(16**l) 

at X=3 %\ -(B-9+3 A %I)/(U8 % I ) 

These values are computed by the residue algorithm given In 
the previous chapter. 

The full algorithm depends on finding the poles of R. 
Mow this Is done will be discussed In detail after consider- 
ing the next Integral. 

K = / R(X) dX 
Jo 

If R(Z) has no real pole which Is positive or zero, 
this integral can be evaluated by Integrating the following 
contour Integral [30] 
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i 



J = I f(Z) dZ, 



f(Z) = PLOG(-Z) R(Z) 



where PLOG denotes the principal branch of the LOG function/ 
I.e., the Imaginary part of LOG(-Z) lies between %?\ and - 
%Pl. The contour C consists of circular arcs of radii a, b, 
and two straight ltnes joining their end-points as shown tn 
fig. 2. 




*-X 



One can verify that contributions from the two circular 
arcs of C vanishes In the limit. Thus, 



= | PLOG(-X) R(X) dX + f PL0G(-X %E ) R(X) dX 

•Jo* 



■fa 



L 



'CO 



at -a» ~ao 

PLOG(X) R(X) dX -/ PLOG(X) R(X) dX - 2 %?\ %\\ R(X> dX 



•/. 



£ 



« - 2 %P\ %\ K 
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Hence the value of the Integral K Is given by 

(3) K = -(£ Res(f) over poles of R Inside .C) 

Therefore/ the value of K Is -1 times the sum of 
residues of LOGC-Z)R(Z) In the complex Z-plane cut along the 
positive real axis. 

It may be of some Interest here to mention that a 
theorem In 1 30 1 gives a value for the integral K which 
differs from the correct formula (3) by a sign. 

Now suppose R(Z) has poles which are real and positive, 
then the Integral K Is divergent. But Its Cauchy principal 
value exists If every such pole Is of order 1. This value 
is given by -1 times the sum of residues of LOG(-Z)RCZ) In 
the entire Z-plane punctured at Z=0, as given In 

(U) (P)K = -( £ Res(f) over poles of R) . 

K diverges if R(Z) has a pole at Z=0. The following 
examples were computed by WANDERER. For ease of reference, 
expression (Dl) Is repeated here 
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2 
X +AX + B 
(Dl) 

k 2 

X + 10 X +9 

(C3) DEFINT(D1,X,0,INF)@ 

2 %P\ B + 6 LOG(3) A + 6 %P\ 
(D3) 

The answer In (D3) Is obtained by locattng the four 
simple poles of (Dl) and applying formula (3). 



(Ci*)l/(X**2+X + 1)@ 
(DU) 



2 
(X + X + 1) 



(C5)DEFINT(D6,X,0,lNF)(a 

2 %PI 
3 SQRT(3) 

Now let us turn to the problem of finding the poles of 
R(Z). The SOLVE [16] program in MACSYMA knows how to solve 
a number of types of equations. This routine Is used In 
obtaining the locations and multiplicities of the zeros of 
Q(x), the denominator of R(X>. The SOLVE program does this 
by factoring Q(X) over the Integers and applying formulas to 
each Irreducible factor of degree less than 5. For factors 
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of higher degree only those of the form aX**n+b will be 
solved. This facility can be augmented by Including more 
special cases for polynomials of higher degree. If the 
integrand Is not rational, the problem of finding poles can 
be much more difficult, although SOLVE can find some of 
these also. In any event, the method of residues can be 
successful only If all relevant poles can be located. 

Let p be a root of Q(X) = of multiplicity m which is 
obtained employing SOLVE. It Is not necessarily true that p 
is a pole of R(Z) of order m. This Is because the rational 
function package [10] in MACSYMA does GCD cancellations only 
over the integers while R(Z) Is a rational function over the 
complex numbers. Therefore p may not be a pole of R(Z) or p 
may be a pole of R(Z) of order less than m. Fortunately 
WANDERER can pretend that p Is really a pole of order m and 
proceed with the algorithm for computing residues which will 
produce a 0, if p is not a pole, and the correct residue, in 
case p Is a pole of order less than m. More discussion on 
the computation of residues can be found in Sect. 6-1. 

Thus, as the reader can easily verify, (D2) of example 
(C2), in page 5? , Is valid even if A=0 and B=l Cor A=0 and 
B=9). As another example, correct answers were obtained in 
the following Integrals even If %\ , the square root of -1, 



6k 



Is treated as a simple pole of the Integrand. 

X - %l 
(D6) 

2 
(X - 2 %I)(X ♦ 1) 

(C7) DEFINT(D6,X,MINF,INF) 

2 %?\ 

(D7) 

3 

(C8) DEFINT(D6,X,0,INF) 

%\ L0G(2) + %?\ 
(D8) 

3 



A special case check Is provided as an auxiliary to the 
above algorithm. 

P-l 
00 X dX 



N M 
° (A X + B) 



P/N /- P/N + M - l\ %PI 

(B/A) 

M - 1 / M 

B N S INC %PI P/N) 



for M 7 N, P positive Integers, M > P/N, P not divisible by N 
and AB > 0. 

For the case M=l, AB<0, the following Cauchy principal 
value Is used. 
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A 1 - P/N %P\ *PI P 

-(- -) (---) COT( ) 

B AN N 

This special case Is not necessary for the algorithm 

but It Is very helpful as far as program efficiency Is 

concerned/ especially for large M or N. 

Sometimes substantial computing time can be saved by 

application of differentiation techniques In Integration 

problems. Consider the Infinite tntegral 

-co 
(5) 



r m dx 

J 2 1( 
"'-* (X + X + K) 



If the algorithm for rational functions were used straight- 
forwardly, differentiation would have to be carried out 18 
times to obtain the sum of residues In case both poles lie 
above the real axis. But since 



dX 
J = 

2 




"* X + X + K 



it Is only necessary to differentiate a quadratic expression 
once to obtain the sum of residues to evaluate the Integral 
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00 dX 



2 
'-* X + X + K 



After this Is obtained/ 9 more differentiations with respect 
to a parameter K are needed. It means a saving of at least 
8 differentiations. 

In fact/ whenever the denominator Q(X) of the Integrand 
Is of the form 

v(X)**n 
with v(X) a polynomial and n s 2, this method of diffe- 
rentiation may be applicable If the degree of Q(X) Is large 
compared to that of the numerator PCX). 

We have Included this technique in WANDERER. Let 

n - deg(PCX)), 
s = deg(v(X)), 
Q(X) = v(X)**m. 

Here Is a brief description of the algorithm. 

If 2 s m*s-n/ this algorithm is not applicable/ 
otherwise (m*s-n > 2) 

(I) Set r to the least Integer 2 (n+2)/s 

(II) If m > r 

(a) Compute the original Integral with m 
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replaced by r and v(X) replaced by v(X) + ZP 
where ZP Is a new parameter Introduced Into the 
problem. Let the answer thus obtained be 
ANS(ZP). 

(b) Return the following as the answer 



m-r 



m-r (r-1)! 



(-1) 



ANS(ZP) 



(m-D! IdZP 



ZP=0 



Cii) otherwise (m=r), the algorithm Is not necessary 
and thus will not be applied. 

The parameter ZP used In this algorithm makes the 
symbolic differentiation In step (ll-b) possible at all 
times. The final answer Is obtained by setting ZP to zero. 
Such a parameter will be used again later and will be 
referred to as the zero parameter. 



For example, 



(D9) 



2 3 
(X + X + 1) 



(C10)DEFINT(D9,X,0,lNF)(a 



k %P\ 



9 SQRT(3) 2 
In this example, we have n=0, m=3 and s=2. Thus r=l 
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and t h e value 




dX 



(X + X + 1) 



I s computed first 
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2. Infinite Integrals Involving Algebraic Irrational Func- 
tions 

2.1 Evaluation by Residues 

Let R(X) be a rational function. Consider 

D = / X**K*R(X) dX, -1< K < 1 

If 

LIMIT(X**(K+1) R(X), X, 0, PLUS) = 0, 

LIMIT(X**(K+1) R(X), X, INF) » 0, 

and R(X) has no poles of order greater than one on the 

positive real axls / then D can be evaluated by applying 
residue theory to the contour Integral 



/. 



f(Z) dZ, f(Z) = (-Z)**K R(Z) 

around the contour C shown In fig. 2. The limit conditions 
given above are convergence tests for D. They are computed 
by DELIMITER. In computing these Integrals/ the Integration 
program will first obtain two sums of residues SI and S2 of 
the function f(Z) by executing the following two steps. 

If R(Z) has poles off the positive real axis, set SI 
to ]TRes(f(Z)) at these poles, otherwise set Si to 
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0. 

If R(Z) has simple positive poles, set S2 to 

ljRes(f(Z)) at these poles, otherwise set S2 to 0. 

Now the answer to the Integral D Is given by 

SPI SI 

%P I COTUPI (K + 1)) S2 

SINUPI (K + 1)) 

For example, MACSYMA produced the following result. 

(CI) l/((l+X)*X**(l/2))@ 

1 
(Dl) 

SQRT(X) (X + 1) 
(C2) DEF1NTU,X,0,INF)@ [3] 

(D2) %P | 

2.2 Integrals Related to the Beta Function 
From the definition of beta function 

BETA(K,S) =/ X**(K-1) (1-X)**(S-1> dX, 



1r 



where K > and S > 0, One can deduce the following rela- 
tion 
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L 



co K " 1 K " S 

00 X D BETAU/S - K) 



dX • 

S K 

(CX+D) C 



w 



here CD > and S > K > 0. 



By a further change of variable In this formula/ the 
following very useful relation Is obtained. 



K - 1 
,a? " BETA(A,B) 



I "_7 



(a) / dX , r > 

S A B 
° (CX +D) C D r 



where A ■ K/r > 0, and B ■ S-A > 0. 

This rather general formula covers many Interesting 
Infinite Integrals Involving algebraic Irrational functions, 
This Is built into the programs by using special purpose 
pattern recognition routines (see chapter 6) to examine 
whether the Integrand Is of the particular form (a). A 
simplification routine for Beta function was also needed to 
put the results In a simpler and better looking form. 



K 
X 



(D3) 



X ♦ 3 

(Ck) DEFINT(D3,X,0,INF)t [12] 

IS THE EXPRESSION 

K + 1 

POSITIVE, NEGATIVE/ OR ZERO 
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POSITIVES 

IS THE EXPRESSION 
- K 
POSITIVE, NEGATIVE, OR ZERO 

POSITIVES 

K 
(Qk) 3 BETA(K + 1, - K) 



2.3 Evaluation by Trigonometric Substitution 

If R(X,Y) Is a rational functions In the two variables 
X and Y, then an Integral of the form 

rao 

U = / R(X, SQRT(X**2-A**2)) dX 
Jo 

can be transformed to an Inftntte Integral of a rational 
function which has been discussed In Sect. U-l. The trans- 
formation below Is easy to verify. 

Let 

2 2 
SQRT(X - A ) 
Y = -. 

X + A 

U becomes 

2 

r* A (Y + 1) 2 A Y Y dY 

i* A / R( , ) 

j 2 2 2 2 

Jo 1-Y 1-Y (1-Y) 
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Now let Z = Y/(l - Y) in the above Integral , It can be shown 
that U Is equivalent to 

2 

f 00 A (2 Z + 2 Z + 1) 2 A Z (Z ♦ 1) Z (Z + 1UZ 

(I) UA/ R( , ) ""' 

/ 2Z + 1 2Z + 1 2 

Jo (2 Z ♦ 1) 

Although the above Is all that Is needed for the given 
Integral U, It Is sometimes more efficient to transform U 
Into a finite Integral of a rational function of the trigo- 
nometric functions. That Is, by setting COS(t) = A/X In U, 
one may deduce that 

r*^ A SlN(t) 

(II) U = A/ RC ,A TANCt)) C )dt 

I COS(t) 2 

Jo COS (t) 

Both of these methods (I) and (If) have been Included 

In the programs. Our heuristic rule for applying these 

transformations Is as follows. 

Apply (II) If It transforms the given Integrand Into 
the form 

C0S(t)**m SlN(t)**n, 
otherwise apply (I). 

For example, 



7k 



1 

(D5) 

2 2 2 
X SQRKX - A > 



(C6)DEFINT(D5,X,0,INF)@ 

1 
(D6) 

2 
A 

2.k Differentiation with Respect to a Parameter 

We have Included two formulas tn WANDERER for Infinite 
Integrals Involving algebraic Irrational functions. They 
are 

"* dX 1 

• — - = U(A,B,C), 

2 3/2 SQRT(C) (B/2+SQRT( AC) ) 
'* (AX + B X +C> 

with A S 0, C > 0, and B > -SQRT(AC), and 

*°° X dX 1 

= _-- = V(A,B,C) 

2 3/2 SQRT(A) (B/2+SQRTC AC) ) 

'" (AX + B X +C) 

where A > 0, C > and B > - SQRT(AC). 

They are Included because any Integral of the form 
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(ill) 



M 

r® x dx 

/ 2 N+ 

J o (AX + B X +C) 



3/2 



can be computed from these known results. Our tool In doing 
this Is differentiation with respect to an appropriate 
parameter. Let Zl, Z2 and Z3 be three zero parameters and 



H = (-1) 



N SQRT(2PI)/2 



GAMMA(- + N) 
2 



The algorithm Is as follows. 



Upon deciding that the Integrand Is as given In (lit), 
WANDERER has obtained the values of M, N, A, B and C. The 
answer to the given integral Is then computed using the 
simple procedure: 

If M B and N = 0, return U, 

If M = 1 and N - 0/ return V, 

If 2N+2 SM, the given Integral Is divergent. 



If N S m, return 
/ 



N 



H< 



d 
N-M 



M 



U(A / B+Z2 / C+Z3) 



dZ3 



dZ2 
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If M ■ N+l, return 



N 

H [—4 V(A,B+Z2,C>, 
l dZ2 / 



If 2(N+1) > M > N+l, then 

(I) If M Is even set r » M/2, return the answer 
/ 



N 



H < 



N-r r 
dZ3 dZl 



U(A+Z1,B,C+Z3) 



(II) If M Is odd, set r = (M-D/2 and return 



N 



H< 



N-r r 
dZ3 dZl 



V(A+Z1,B,C+Z3)> 



/. 



A s Imple example Is 



(D7) 



(X + X + 1) 



15 
2 



<C8) DEFINT(D7,X,0,INF)@ 

1021* 
6567561 



(D8) 
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3. Infinite Integrals Involving Trigonometric Functions. 

3.1 Integrals from Minus Infinity to Infinity 

Let R(X) be a rational function which has no real poles 
and 

LIMIT(R(X),X,INF) » 0, 

The Integral 



J-to 



°° i%\ m X) 



L « / %E R(X) dX, m > 

J-tD 

Is convergent. Its value can be obtained by evaluating the 
contour Integral 



\%\ m Z) 
J = / %E R(Z) dZ 



'1/ 



around the contour C given In fig. 1. 

Let Cr be a circular arc with center at Z = 0, radius r 
and argument t, tl < t < t2. Jordan's Lemma [7] shows that 
If f(Z) approaches uniformly on Cr as r approaches INF, 



f (%l m Z) 
IMIT %E 



LIMIT I |E f(Z) dZ»0, for tl 3 0, t2 S ?PI 



and that 
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IMIT / 



(-%\ m Z) 
LIMIT/ %E f(Z) dZ=0, for tl2£$;PI, t2 S 2%PI . 

Therefore L « J as r tends to Infinity and the method 
of residues can be applied. Moreover, for m < one may 
use the same method by using a contour similar to C In the 
lower half complex Z-plane. 

For real -m, n and p, let T(X) be SIN(mX), C0S(nX> / 
%E**(.%\ p X) or a function Involving sums and/or products of 
these functions. In complex exponential form T(X; Is a sum 
of constant multiples of functions of the form %E**UI k X), 
k real. 
Therefore 

I T(X) R(X) dX 

can be Integrated using complex contour Integral. In doing 

such a problem/ terms In the Integrand are sorted into two 
parts. One part requires a contour In the upper-half Z- 

plane, the other a contour In the lower-half. 

Before discussion of additional methods, let us see a 
few examples computed using methods discussed so far. 
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(Dl) 



(C2) DEFINT(D1,X,MINF,INF)Q 



(D2) 



COS(X) 



X ♦ 1 






[12] 



(C3) X*SIN(X)/(X**2+1)@ 
(D3) 

(C4) DEFINT(D3,X,MINF,INF>Q 

(C5)X*COS(X)/(X**2+l)@ 
(D5) 

(C6) DEFINT(D5,X,MINF,INF)@ 

(D6) 

(C7) X*COS(X)/(X**2+X+l)@ 

(D7) 

(C8> DEFINT(D7,X,MlNF,INF)@ 



X SIN(X) 

2 
X ♦ 1 



%PI 

%E 

X COS(X) 

2 

X + 1 



X COS(X) 



X + X + 1 



[12] 
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(D8) 

1 (SQRT(3) - %|>/2 

WI/2) (( %\) %E 

SQRTC3) 

1 (- SQRK3)- %\)/2 
+ (II ) %E ) 

SQRTC3) 

(C9) 1/UE**(*!*X)*(X**2 + 1))@ 

1 
(D9) 

II X 2 
%E (X + 1) 

(CIO) DEFINT <D9,X,MINF,INF)@ 

%PI 

(D10) 

%E 
(Cll) SIN(X)*D9@ 

SIN(X) 
(DID 

%\ X 2 
*E (X + 1) 

(C12) DEFINT(D11,X,MINF,INF)P 

%?\ 
(D12) 

*E 

(C13) D11*(%E**(-%I *X>) 

SIN(X) 
(D13) 

2 %\ X 2 

*E (X + 1) 

(C1U) DEFINT(D13,X,MINF,INF)@ 
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(D1U) 

2 



3.2 Integrals from to Infinity 

Let us now discuss some Integrals Involving trigonome- 
tric functions over the range (0 INF). First consider the 
Integral s 



■/.' 



n 
ll=| C0S(k X ) dX 

o 



(-co n 

and 12 A SINCk X 



CO 

) dX 

o 



where k a nonzero real constant and n > 1. 

Here again/ the residue theory can be applied to 
evaluate the Integrals. This time the shape of the contour 
Is slightly different. The contour Is a sector of a circle 
consisting of a portion of a circular arc with center at the 
origin and two straight lines joining the end points to the 
origin, (fig. 3) The sector angle depends on n and Is 
%PI/(2n). 

Without loss of general I ty, let us assume that k > 0. 



83 



Therefore, taking real and Imaginary parts of the above 
relation, we have 

II = C0SUPI/2n) G 
and 12 = SINUPI/2n> G. 



For example 

7/3 
(D15) C0SC9 X ) 

(C16) DEFINT(D15,X,0,INF)@ 



3 3 %PI 

3 GAMMA(-) COS( ) 

7 Ik 

(D16) 

3/7 
7*9 



(C17) SIN(9*X**(7/3))@ 

7/3 
(D17) SINC9 X ) 

(C18) DEFINT(D17,X,0,INF)@ 

3 3 %P\ 

3 GAMMA(-) SINC > 

7 lh 

(D18) 

3/7 
7*9 

Incidentally, if n as 1, II and 12 do not converge, 
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A more general Integral than II and 12 is 



Jo 



m n 

X EXPUI k X ) dX 



where n > 0, Rl (m) > -1/ k real and nonzero and n-Rl(m) > 1. 

We shall consider the case k < 0. The case k > Is 
entirely analogous. Let us take a sector-shaped contour In 
the fourth quadrant as shown In fig. k. 



<** 




'».4» 



-# 



By Cauchy Integral theorem/ 



I 



n 



w 



hlch Impl les that 



m Ul k Z ) 
Z %E dZ » 



-co 



m n 

J = EXPC-%PI %l(m*l>/(2n>>/ R EXP(k R ) dR 



85 



s 
* EXP(-*PI %l(m+l)/(2n)) GAMMA ( ) 

s 
n(-k) 



where s - (m+l)/n. 



Let us look at an example. 



(2 %\ + 3) 
X 
(D19) ■ 

3 

(21 X/2) 
*E 

(C20) DEFINT(D19,X,0,INF)<a 



2 %\ + k 

GAMMA( ) 

3 
(D20) 

(2 %\ *k) I 3 
3 (21 / 2) 

A similar result for k > Is given by 



(m+l)/n 
J = EXP(%PI %l(m+l>/(2n)) GAMMA( ). 

(m+l)/n 
k n 

The following Integrals can be obtained readily from 

the above relation. 



n m r n m 



COS(k X ) X dX and / SIN(k X ) X dX 
o 

For Instance, 



Ir 
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■» r -s r (r - 1) f* r-2 -Cs-2) 
SIN (X) X dX = / SIN (X) X dX 

(s-l)(s-2)Jo 



r f 00 r -Cs-2) 
+ SIN (X) X dX, 

(s-l)(s-2)Jo 



where r > Cs-1) > 1. As one can see / repeated application 
of this relation will reduce U to a sum of Integrals of the 
same form as U but with N = 1 or 2. For N=l, AUDITOR uses 
the following formula 

p -1 / P-1 
SIN (X) X dX =1 SIN (X) dX, P > odd. 

O Jo 

The Integral on the right hand side can be evaluated easily. 
Integrals of this type are considered In Sect. 5-2.2. For 
N=2, we have 



q -2 %PI /q-3/2\ 
SIN (X) X dX = — 

2 I q-1 




w 



here q is an Integer 2 2, For example 



SIN(R X) 

(Dl) 

X 



(C2) DEF INKS ,/;,(;, IMF) a 
(D2) 



%? 



( C 3 ) S I N ' ■' ■ * X ) * * 2 / ( X * * 2 ) @ 



(D3) 



SIN (Q X) 
2 



(ClO DEFINT($,X,0, INF)@ 
CD4) 



2 
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k. Infinite Integral Involving Logarithm Functions 

4.1 Evaluation by Recursion and Contour Integration 

Let R(X) be a rational function which Is even In X. 
The Integral 



-co 



N 



V(N) =/ LOG (X) R(X) dX, R(X) = P(X)/Q(X), 

with deg(Q(X))-deg(P(X) ) 2 2, can be evaluated by applying 
residue theory In a recursive manner. Recall that LOG Is an 
abbreviation of PLOG when the argument Is real and positive. 
Consider the contour Integral 



JC 



N) « / 

k 



N 
F(Z) dZ, FCZ) = PLOG (Z) R(Z) 



where C Is the Indented contour In fig. 5. As r approaches 




► X 



f^JL 
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1) Declare each of AV and AJ to be an 1 X N array 

2) Set I to then go to step k 

3) Set AJ(i) to J(I)/2 

k) If I = N, return the answer V(N) 

5) Set AV(!) to V( I) 

6) Set ! to 1+1 then go to step 3 



Computation of V( I ): 



If I = 0, return the value of the Integral 

-co 

R(X) dX, 



(storing away Information about poles and 
corresponding residues of R(Z) for computation of J), 
otherwise, compute from the formula 



V(l) » AJ(t) - - 
2 




K K 
%P\ %\ AV(I - K), 



Computation of J(I)/2 : 

Using existing Information of poles and residues of 
R(Z) compute by residue theory from the formula 

%P\ %\{ £ Res(F) over poles of R Inside C). 
For example. 



92 



(Dl) 



LOG (X) 



X + 1 



(C2) DEFINT(D1,X,0, I N F ) @ 



[11 



(D2) 



8 



k.2 Method of Differentiating and Introducing Parameters 

A very useful method for Integrands Involving LOG(X) Is 
differentiation. Consider an Integral of the form 



-co 



K 



A(K) = R(X) X LOG(X) dX, K t and -1 < K < 1, 

where R(X) Is rational In X and k Is a parameter which 
occurs nowhere else in the Integrand. If the Integral 



r cD K 

B(K) = R(X) X dX 



Is convergent and can be evaluated, then A(K) Is given by 



A(K) = -- B(K). 
dk 

This method of differentiation is valid because B(K) is 

convergent and A(K) Is uniformly convergent for every closed 

interval contained In the set of points (K | MO and -i<K-=l). 



For example 
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(D3) 



(CO DEFINT(D3,X,0,INF)@ 



X LOG(X) 
X + 3 



112] 



K K D 

(DIO L0G(3) 3 BETACK+1,- K) + 3 (~BETA(r, - K) 

Dr 
D 

BETAU + l,r) 

Dr 



r»(l+K) 

) 

r=-K 



In the foregoing discussion, K has been assumed 
symbolic. Yet, In solving actual problems K may very well 
be a number. This difficulty can be overcome by Introducing 
a zero parameter ZP. We first replace X**k by X**(k+ZP) In 
f(X). This permits us to proceed as above and then diffe- 
rentiate with respect to ZP. After differentiation the 
result Is then evaluated at ZP » 0. 

The property of the logarithm functton 

-LOG(X) = L0GC1/X) 
can sometimes be used In evaluating Integrals of the form 



■r 



LOG(gCX)) f(X) dX 
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Is convergent, then 



J = L0G(g(X)) U(X) 



U(X)g'(X> dX 
( ). 

a j a g(X) 



The SIN Integration program [20] In MACSYMA Is used to 
compute the Indefinite Integral U(X). WANDERER uses this 
method to evaluate Integrals of the form 



-<"-L a 
X LOGU + X ) dX, 



a + 1 > L > 1. 



Here Is an Integral evaluated by the method of Integra- 
tion by parts. 



(D7) 



7/2 

LOGU + 1) 



(C8) DEFINT(D7,X,0,INF)Q 

*PI 



(08) 



5 %P\ 
SINC ) 

7 
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5. Infinite Integrals Involving Exponential Functions 

5.1 Method of Substitution 

Integrals with Integrands which are rational functions 
of £E**(K*X), K real and non-zero, are relatively easy to 
compute. Without loss of generality/ let us assume K > 
and consider 

(1) / RUE**(K*X)> dX. 

J-co 

If one makes a change of variable 

(2) Y » *E**(K*X), 
then (1) becomes 

GO 

(1/K) I R(Y)/Y dY. 



f 

r K> I R< 



This Integral converges If (1) does. That Is to say 
R(Y) has Y as a factor If (1) converges. For the Integral 

-co 

(3) / RUE**(K*X)) dX 



r 



one can make a slmtlar change of variable 
U) S + 1 = %E**(K*X) 

which gives 
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r 

(5) / R(S+1)/(S+1) dS. 

Integrals of the form of (3) and (5) can be evaluated 
by contour Integral and residue theory as Indicated In sect. 
1. The methods (2) and (k) can In general be applied to any 
Integral of the form 

/ f(g(X)) dX 

and the resulting Integral may be much simpler. For 
Instance 

1 
(Dl) 

X/3 5 
*E ( + 7) 

X/3 

2E 

(C2) DEFINT(D1,X,0,INF)Q 

7 
3 LOGC — ) 

12 
(D2) 

5 

XA 
%E 
(D3) 

X/2 

9 ?E + k 

(C«0 DEFINT(D3 / X / M!NF / INF)@ 
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(DlO 

3 



5.2 Use of Contour Integration 



A somewhat more Interesting Integral Is 

CO 



= / R{% 
J-cw 



E**X) PCX) dX 



where P(X) Is a polynomial and R(X) a rational function with 

complex coefficients such that 

UM!T(R(%E**X),X,INF> = 
and 

UMIT(RUE**X),X,MINF) « 0. 

Let us first determine a polynomial Q(Z) with complex 
coefficients such that 
(6) Q(X) - Q(X + 2 %?\ %\) = P(X). 

Q(Z) exists and can be computed by the method of 
undetermined coefficients. Now consider a contour Integral 



i 



J = / RUE**Z) QCZ) dZ 
6 

taken around a rectangular contour as shown In fig. 6 
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VK 



+-X 



As the absolute value of X approaches Infinity, 
contributions from the vertical segments of the rectangular 
contour vanish. Hence, 



CO 

J = / R(^E**X) Q(X) dX 



'-co 



-/ R(% 
J-co 



E**X) Q(X+2 %?\ %l) dX. 



It follows from (6) that l=J. Mow the problem of evaluating 
the Integral I has been reduced to finding the residues of 
R(2E**Z)Q(Z) for se lm(Z) < 2 %Pt . To do this, the poles 
of R(Z) are obtained first. If to Is such a pole then 

GLOG(w) 
where GLOG stands for the branch of LOG with Imaginary part 
between and 2*%PI, Is a pole of the same order for 
R(%E**Z)*Q(Z) Inside the closed contour. All such poles can 
be obtained In this manner. One example Is 



(D5) 



SINH(X) - %\ 



(C6) DEFINT(D5,X,MINF,INF)@ 
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(D6) %PI 

Note that the Integrand (D5) has a pole of order two at 
Z = %PI %\/2. 

5.3 The GAMMA Function and Related Integrals 

A very Important function closely related to the 
evaluation of Infinite Integrals Involving exponentials Is 
the GAMMA function generally defined as 



■F. 



-co 

GAMMA(Z) =| 2E**(-t) t**(Z-l) dt, R1(Z) > 0. 

'o 



Also of use Is Its logarithmic derivative/ the PS I 

function 

d 
PSKZ) »(— GAMMA(Z) ) /GAMMA (Z). 
dZ 

A simplification routine for GAMMA function has been 
written to make use of the many properties of this function, 

From the definition of the GAMMA function/ one can 
derive the following very useful relation 



-00 D Q 

D C - A X GAMMA(a) %E 
X 2F. dX » » 

a 
A B 
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w 



here a = (D + 1)/B, R1(A) > 0, R1(D) > -1 and R1(B) 



nonzero. 

WANDERER has programs designed to recognize this form 
and return the result. Of course these programs have to 
examine the signs of the relevant quantities carefully 
before generating an answer. 

For example: 

2 
-X 

(D7) %E 

(C8) DEFINT(D7,X,0,INF)@ 

SQRT(fcPI) 



(D8) ■ 

2 



5.k Integral Related to the Laplace Transform 

Let f(t) be a function of a real variable t, then Its 
Laplace transform L(f(t)) Is defined as 

"• C-t s) 
L(f(t)) » F(s) =| %E fCt) dt. 



f ( " t 
Jo 



Many of such Integrals can be evaluated by the pro- 
grams. Here are some examples produced by WANDERER. 
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CHAPTER V 
INTEGRALS OVER A FINITE RANGE 

0. Introduction 

Many finite Integrals are proper Integrals whose 
Indefinite Integrals exist In closed form and can be com- 
puted rather easily. For such an Integral/ the evaluation 
method is very straight-forward. WANDERER simply obtains 
the corresponding Indefinite Integrals and then substitutes 
the limits of Integration. The ant Iderlvat Ives are computed 
by use of SIN. This method will be referred to as the antl- 
derlvatlve method. 



A finite Integral 



r 



b 
f(X) dX 

a 



Is Improper If f(X) becomes Infinite at some point c, 
b s c & a. In order to avoid having to spend computation 
time looking for an unknown number of singularities of f 
between a and b, WANDERER will assume that Integrands of 
finite Integrals, with the exception of rational Integrands/ 
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can become Infinite only at the end points of the range of 
Integration. This Is not a severe restriction and does not 
decrease the number of Integrals It can handle. This Is 
true In the sense that any given range can be subdivided to 
conform to the above convention. This convention makes It 
easy to determine whether a given finite Integral Is 
Improper. WANDERER simply checks the value of the Integrand 
at the limits of Integration. If the given Integral Is 
Improper, Its divergence Is tested before any attempt at 
evaluation. WANDERER uses a limit test for absolute 
divergence which Is discussed In Scet. 6-6. If the given 
Integral diverges, WANDERER will so Indicate In the output. 
If the antlderlvatlve can be computed, then the answer Is 
sometimes obtained by employing a limiting process when 
substituting the upper and lower limits of Integration. 

In this chapter, attention will be focused on the 
definite Integrals whose corresponding Indefinite Integrals 
are difficult to compute or do not exist. 

1. Finite Integrals of Rational Functions. 
For Integrals such as 



I r(x: 



U = R(X) dX 

'a 
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where R(X) Is rational In X, a and b are flnlte / WANDERER 

computes U by transforming It Into an Infinite Integral by a 
change of variable. 

Let us write 

(X - a) (a + Y) 

(1) Y = -- , X / 

(b - X) (1 ♦ Y) 



(b - a) dY 
dX = 

2 
(1 + Y) 



then It Is evident that 

-co 



I b Y + a 
RC , 

J« Y ♦ 1 



dY 
U = (b - a) / R( , ). 

(Y ♦ 1) 



This Integral can be Integrated readily by methods of con- 
tour Integration and other means discussed In sect. k-1. 

As an example let R(X) be the expression (Dl) 

1 

(Dl) 

2 
X - 3 

The Indefinite Integral of R(X) from to 1 computed by 

the substitution method above Is given by 

(C2) DEFINT(D1,X, 0,1)8 
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LOG(2 - SQRT(3)> 



(D2) 

2 SQRT(3) 

One knows that the Indefinite Integral of (Dl) exists. 
In fact the following has been obtained using the command 
INTEGRATE In MACSYMA 
(C3)INTEGRATE(D1,X)@ 

(D3) (LOG(2 X - 2 SQRT(3)> - L0G(2 X + 2 SQRT(3>)) 

/(2 SQRT(3)) 
The reader nay easily obtain a result equivalent to 
(D2) by Substituting In (D3) the limits of Integration. 

2. Rational Functions of Trigonometric Functions 

2.1 A Typical Application of Contour Integration 

If R(X,Y) Is a rational function In two variables X and 
Y, an Integral In the form 



U = / R(COS(X),SIN(X)) dX 
■Jo 



is easily transformed to an Integral around a closed con- 
tour. By setting 

2 2 

Z + 1 Z - 1 

COS(X) = , SIN(X) / 

2 Z 2 *l Z 
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dZ 
and dX a ----, 

%\ Z 

the given Integral U becomes 

2 2 



if Z+1Z-1 dZ 
- I F( , ) — 

%\JM 2 Z 2 %\ Z Z 



where CI Is the positively orientated unit circle with 
center at Z = 0. This contour tntegral can then be 
evaluated by finding the sum of residues Inside the circte. 

Actually this transformation can be applied In general 
to any Integral In the form 



L 



R(%E**(%I X)) dX 



where R Is rational, by the change of variable 

Z = 2E**(%I X) 

The transformation process Is simple. The key point In 
this algorithm Is the determination of whether the Integrand 
Is In fact a rational function of %E**(%l X). 

In calculating the sum of the residues, only poles 
Inside the contour contribute. Simple poles on the unit 
circle, I.e. those with absolute value 1, cause the 
principal value to be computed. The Integral is divergent 
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If poles of order greater than 1 happen to be on the circle. 

Examples: 

iCk) C0S(X)**2-SIN(X>@ 

2 
(Dlf) COS (X) - SIN(X) 

(C5) DEFINTCDt,X,0,2*3PI)@ 

(D5) %PI 

(C6) %E**(2*%l*X)/(%E**C%l*X)+3)@ 

2 %\ X 

(D6) 

%\ X 
%E + 3 

(C7) DEFlNTCDe.X.O^nPDQ 

(D7) 2 2PI 

2.2 Utilization of the Periodicity of the Trigonometric 
Functions 

In this section Integrals of functions Involving trigono- 
metric functions over a variety of ranges will be 
cons tdered. 

Let T be a function of X defined by 



M N 
T(X) = COS (X) SIN (X) 
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where R1(M) > -1 and Rl(N) > -1. The following formula can 

be deduced from the definition of Beta function. (Sect, k- 

2.2) 

M + 1 N + 1 
^ l9 GAMMA C ) GAMMA< ) 

"*'* 2 2 



(1) / TCX) dX » 

2 GAMMAUN + M + 2)/2) 



A simple example Is 

1/3 1/2 
(D8) COS(X) SIN (X) 



(C9) DEFlNT(D8,X,0,%PI/2)@ 

2 3 
6 GAMMA(-) GAMMA(-) 

3 h 
(D9) 

5 
5 GAMMA (--) 
12 



The usefulness of (1) Is Increased by the fact that It 
Is possible to express definite Integrals of T(X) over a 
variety of ranges In terms of that of TCX) over (0 2PI/2). 
For Instance, the following relations are true for any func- 
tion f. 

f(SlN(X),COS(X)) dX =» 



L 



Ti/Z 
(f(SI.N(X),COS(X)> + f(SIN(X),-COS(X))) dX . 



Ill 



f(SlN(X),COS(X)) dX = 

f(SIN(X),COS(X>> dX+ / f(-SIN(X),-COS(X)) dX. 
O Jo 

More general ly, let f(X) be a periodic function of X 

with period 2%P I . That Is 

f(X + 2%P I ) = f(X). 

An integral of f(X) over some range (a b) / 

-b 



r 



s 

'a 



f f(x) dx, 
J a 

can always be written as a sum of Integrals In the form 

(2) n f f(X) dX + f f(X) dX - f f(X) dX 
Jo -Jo -*o 

for some Integer n and 1%?\ > c, d S 0. This Is true for 
there exist integers p and q such that 

a = 2 p %?\ + d, 
and b ■ 2 q %?\ + c. 

Then S is equivalent to the sum (2) with n = (q - p). 

Programs have been written to perform this reduction 
and they are applied when the Integrand has a period 2%?\ 
and the difference (a - b) has %?\ as a factor, some 
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examples computed by WANDERER are 

1/3 2 
(D10) COS(X) SIN (X) 



(Cll) DEFINT(D10,X,-*PI/2,!&PI/2)@ 

2 

18 SQRT(^PI) GAMMA(-) 

3 
(DID -—• - 

1 
7 GAMMA(-) 
6 

(C12) C0S(X)**3*SIN(X)**2@ 

3 2 
(D12) COS (X) SIN (X) 

(C13) DEF INK D12,X,3*3;P I /2,3*%P I) (3 

2 
(D13) 

15 



3. Finite Integrals of Algebraic Irrational Functions 

3.1 Rationalizing the Integrand 

If R(X,Y) Is a rational function In X and Y, the Inte- 
gral 



2 2 
(1) K = I R(SQRT(A - X ),X) dX, 



i- 



where A ^c > b S -A, can be rationalized. 
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Let us write 

2 2 
A - SQRT(A - X ) 

(2) Y = — 3 v(X), 

X 

which gives 

2 A Y 

(3) X = 

2 
1 + Y 

Substituting (3) for X, (1) becomes 



2 2 

VCC) A(l - Y ) 2 A Y Cl-Y ) dY 



J . 2 2 2 2 

J Vt*> 1 + Y 1 + Y (1+Y) 

which Is an Integral of a rational function. Integrals of 

this type have been discussed In Sect. 1. 

Similarly for the Integral 

f C 2 2 

J =* / R(X,SQRT(X - A )) dX 

Jb 

where c > b ^ A > 0/ the change of variable 

2 2 
SQRKX - A ) 

t = - = u(X) 

X + A 

can be made to convert J to 



1 



lit* 

uxc) 2 

A Ct + 1) 2 A t t dt 



k A / R( , ) ■ 

UCM 2 2 2 2 

UCW 1 - t 1 - t (1 - t ) 



which Is also an Integral of a rational function. Here are 

two such Integrals evaluated by WANDERER. 

1 
(Dl) 

2 
X SQRT(X - 9) 

(C2> DEFINT(D1,X,3,«0@ [12] 

-%\ 3 + %\ SQRTC7) 
(D2) — PL0G( r) 

3 k 

(C3) l/((X+l)*U-X**2)**(l/2))@ 

1 
(D3) 

2 
(X + 1) SQRTU - X ) 

(CU) DEFINT(D3,X,0,2>@ 

L0G(2 + SQRTC3)) 
(Dl») 

SQRT(3> 

Another method Is to try to transform the given Inte- 
gral to an Infinite integral by the change of variable given 

In (1) of Sect. 1. For Instance the Integral 

•B 

RCX,SQRT((X - AXB - X)>) dX 
'A 
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can be converted to 



-co 



(B + A) Y (B - A) Y dLY 

(B - A)/ R( / ) 

Y + 1 Y ♦ 1 2 

(Y + 1) 



The method of rationalization can also be applied to 



I 



R(X, (CX + D)**(l/Q)) dX 



where Q Is an Integer / C and D are constants, and the range 
of integration needs not be ftnlte. 

For Integrals of this type the substitution 

1/Q 
Y = (C X + D) , 



Q 
Y - D 
X = / 

C 

Q - 1 
Q Y dY 
dX / 

C 



will convert the given Integral Into that of a rational 
function In Y. 
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The principal task In this conversion method Is to 
recognize that the given Integrand Is of the particular form 
and to obtain the Integer Q. in essence, the algorithm for 
doing this Is 

(I) Obtain a list L of all distinct Irrational 
parts In the given Integrand. 

(ID If elements In L are fractional powers of an 
Identical linear polynomial In X, the pattern Is 
matched and Q Is set to the 1cm of the denominators 
of all the exponents of the linear polynomial, 
otherwise the pattern Is not matched. 

3.2 Integrals Related to the BETA Function 

The BETA function Is defined by the Integral 



f l K - 1 



L - 1 
BETA(K,L) = / X (1 - X ) dX 



for RKK) > and Rl(R) > 0. From this definition one may 
readily deduce the relation 



f V - 1 

■Jo 



C L - 1 1 K 
(3) / X (1 - X ) dX = - BETA(-,L). 

• o C C 

WANDERER applies this formula by recognizing the form 

of the given Integrand. In chapter 6 some techniques and 
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programs for pattern recognition are discussed. For Inte- 
gral s In the form 



K - 1 C C L - 1 
(X - A) (B - X ) dX, 



a simply substitution, 



Y = (X - A)/B, 
will transform it into (3). 
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I*. Finite Integrals Involving Logarithm Functions 

When the Integrand Is a function of LOG(X), a given 
Integral may/ In many cases, be evaluated by transforming It 
into one which Involves exponential functions. Consider 

b 
f(L0G(X)) dX, b > a 5 0, 




the substitution 



(I) 


-Y 

X » %E 


or 
(II) 


Y 
X » %E 


converts U to 





for 1 > a, 



for a 3 1- 



/. 



* -Y 

f(-Y) %E dY, r= -LOG(b), s= -LOGCa), 



T 
or 



l 



s Y 

f(Y) |E dY, r=L0G(a), s= LOG(b), 



'T 

respectively. In case a = and b = 1 (or INF), the use of 
(I) or (II) will result In an Infinite Integral which can 
often be evaluated readily using methods provided In 
WANDERER. For Instance by use of (I) the Integral 



I 



L0G(X ) dX 
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becomes 

-co 



L 



q "Y 
(-Y) *E dY, 
o 



which Is an Infinite Integral already studied (Sect. U-5.3). 
The method (I) or (II) can also be applied to Integrals of 
the form 



r 



r 
X f(LOG(X)) dX, 

t ' a 

as In 

SI R f 09 S -RY-1 

LOG (-) X dX = / Y %E dY. 

X 



f« 



For Integrands involving LOG(f(X)), It may sometimes be 
possible to simplify the given Integral by the substitution 

Y = f(X), X » g(Y) 
where the function g Is the Inverse function of f. A very 
simple application of this method Is conversion of 






CO 

X 2 L0G(SQRT(X)+a) dX 



to the Integral 



f CO 



2 / (Y-a) LOG(Y) dY. 



Techniques of differentiation with respect to a 
parameter, as detailed In Sect, k-k.2, can be employed for 
finite Integrals with a factor 
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K 
X LOGCX), RKK) f 0, 



In the integrand. 

Examples : 

K 
(Dl) LOG (X) 



(C2) DEFINT(D1,X,0,1>@ C 12 l 

K 
(D2) (-1) GAMMAU + 1) 



(C3) LOG(X)**(l/2)/(X**2)@ 

SQRT(LOG(X)) 
(D3) 

2 
X 

(CIO DEFINT(D3,X,1,INF)@ [12] 

SQRTUPI) 

(D»0 

2 

(C5) X**(l/3)*(-LOG(X))**(-l/2)@ 

1/3 
X 
(D5) 

SQRT(-l) SQRT(L0G(X)) 

(C6) DEFINT(D5,X,0,1>@ 

SQRTUPI) 
(D6) 

k 
SQRT(-) 

3 
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(C7) LOGCX}*( ( i-/.**i l/2))/X)**CI/2)Q 

UQRTQ - SQRT(X)) LOG(X) 

(D7) --~ 

SQRT(X) 

(C3) DEFifiKD/,, a, - : . t I) 4 



(D8) 



D 3 | 

(~ BETACr,-)) 
Dr 2 



r = l 



123 



(I) Set V to the quotient of V(Z)/(Z - p)**m 
which is computed by long division, 
(ii) Return the result computed from 



(m - 1)! 



m-1 
/dj U(Z) 
dZ 



Z = P 



V(Z) 

otherwise/ (m > 1, V not a polynomial) compute and 
return as answer the limit 

m-1 



L 
Z 



-p (m - 1)1 IdZJ 



m 



(Z - p) f(Z) 



1.2 Evaluation of Contour Integrals by Residue Theory 

In order to evaluate an integral of f(Z) around a 
closed contour C by residue theory/ it is necessary to 
locate all poles of f Inside C. After this Is done, the 
remalntng problem Is to compute the sum of residues of f at 
these points efficiently. The difficult part Is finding 
poles. WANDERER employs the SOLVE routine In MACSYMA to 
solve V(Z)=0. 

For V(Z) a polynomial In 1, SOLVE finds its zeros by 
factoring over the tntegers and applying formulas to each 
factor of degree less than 5. For factors of higher degree 
only those of the form a*Z**n+b will be solved. The problem 
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of factoring large arbitrary polynomials Is non-trlvlal, to 
say the least. The development of a more powerful factoring 
algorithm which will factor polynomials over a larger ring 
than the Integers would certainly be helpful to SOLVE and 
WANDERER. The location of zeros Is usually more difficult If 
V(Z) Is not a polynomial In Z but some more complicated 
function. In such a case/ WANDERER usually uses methods 
other than the residue theory. An exception Is when V(Z) is 
a polynomial In %E**Z which has been discussed In Sect, k- 
5.2. 

Thus within the limitations of SOLVE, V(Z) will be 
solved and Its zeros sorted Into a list of pairs. Each pair 
containing a zero and Its multiplicity/ such as 

L = ( (Zl/ml) / (Z2/m2) / (Z3/m3) / . . . ). 

L Is then sorted Into two lists LI and L2/ discarding 
poles outside the closed contour C/ such that 

LI » a list of all simple poles 

L2 = a list of other poles paired with their 

mul tlpl tcl ties. 

At this point/ we can apply the RESIDUE ALGORITHM to 
obtain the desired sum. To avoid repeated calculation of 
V'(Z), the program checks whether LI Is empty. If LI Is not 
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an empty list, V'(Z) will be computed and stored for 
possibly repeated reference later In the computation. 

Note that it has been assumed from the beginning of 
this section that a pole of f(Z) would be a zero of V(Z). 
This assumption Is quite reasonable for almost all of the 
applications In evaluation of definite Integrals of 
elementary functions. However a zero of V(Z) need not be a 
pole of f(Z>. For Instance Z - %\ Is a zero of 



Z - 1 



but not a pole of 



(Z + %\) 

k 
(Z - 1) 



For each zero p of V(Z) we may check the value of U(p) 
to see If p Is really a pole of f(z). Although It Is not 
clear what can be done If U(p) = 0, since p may still be a 
pole of lower order. A better method Is to ignore the fact 
that If V(p) = 0/ U(p) may also be and pretend that p Is 
an actual pole of f(Z). This Is valid because a residue at 
any removable singular point will turn out to be 0. 
Furthermore, the residue obtained at a pole of order m Is 
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2. Obtaining Real and Imaginary Parts 

In the course of evaluating a definite Integral by com- 
plex contour Integration, the need to take the real or the 
imaginary part of an expression often arise. For Instance, 
to see whether a pole, p, lies above or below the real axis 
the sign of Im(p), the Imaginary part of p, Is examined. To 
determine If p lies Inside the unit circle at Z=0, It Is 
needed to compute ABS(p) which Involves taking the real and 
Imaginary parts of p. 

The algorithm for obtaining Rl(p) Is presented as a 
representative of similar procedures used. 

Algorithm REALPART(p) : 

1) If p Is %\ return 0, If p Is a number or any 
other atomic symbol, return p. 

2) If p Is a sum, (p » £ p j ) then return 

T REALPART(pj) 

3) If p Is a product (p a pl*p2) then return 

REALPARTC pl)*REALPART( p2 > 
- IMPART(pl)*lMPART(p2> 
i*) If p = ?E**pl return 

?E**REALPART(pl)*COS(IMPART(pD) 
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5) If p = pl**p2, convert p to 2E**(p2*LOG(pl) ) then 
go to step i». 

6) If p = LOG(pl), return L0G(ABS(pl)) 

7) Otherwise, return the form Rl(p). 



3. A Heuristic Pat fern Recognition Progr 



am 



It Is often the case that some pattern recognition Is 
needed, at one stage or another during the evaluation of an 
Integral. Although many tntegrals can be evaluated without 
any pattern recognition, this capability remains Important 
to WANDERER. One specific pattern shall be discussed as a 
representative of such methods In WANDERER. Consider the 
pattern 

N M 
PCX) « (B X + A) 

where B, N, A and M are free of X (I.e. do not Involve X) 

and all except A must be non-zero. This pattern should 

match, for example, every one of the expressions 

1/2 -1/3 
X, 2 X + 1, (X + %\) 



P q 2 

(X - 1) , and X +2 X + 1^. 

The last expression Is the expansion of 
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2 

(X + 1) . 

Expressions In expanded form present some trouble for 
the recognition algorithm. This difficulty Is overcome by 
the use of differentiation and rational simplification. 
That Is If an expression E(X) Is equivalent to an expression 
matching the pattern PCX), we can compute 

E' (X) / E(X) = U(X) / V(X), 
cancelling all common factors In the numerator U(X) and 
denominator V(X). Then V(X) should match p(X) with M =» 1. 
As a result of this match, some values are assigned to the 
variables A and B. By use of these values the correct value 
of M and N can be recovered from U(X). The values for A and 
B thus obtained may differ from the true values by a 
constant factor. This happens whenever these true values 
have common factors. The real values of A and B can be 
determined by comparing A**M to E(0). This procedure was 
suggested by Moses. 

The full algorithm used for matching P(X) will be 
described. Let It be called PM . PM uses another routine PfJ 
which recognizes the pattern 

N 

PCX) - (B X + A) 

As one can see in the following algorithm, the values 
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set 



1/m 1/m 

A - A r and B n B r, 

then P(X) Is matched. 



Algorithm PN(ECX),X) 

1) |f E Is free of X, the pattern not matched. 

2) Match A to ECO) and set E to E - A. 

3) If E Is of the form r*X**s, the pattern Is 
matched (B to r, N to s). Otherwise there Is no 
match. 
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k. A Procedure for Change of Variables 

Substitution of a new variable Y for a subexpression/ 
say g(X), of the Integrand In a given Integral Is a 
frequently used method In Integration. Let the given Inte- 
gral be 



■r 



f(X) dX, 
a 

then the transformed Integral would be In the form 



■f 



J =1 F(Y) dY 
'c 



where F(Y) might be considerably simpler than f(X). FCY), c 
and d are computed by a procedure which Is called whenever a 
change of variable Is needed. It makes use of two other 
modules of MACSYMA, namely SOLVE and DELIMITER. 

I) Use SOLVE to solve for X In Y = g(X), obtaining 

X = h(Y), the Inverse function of g. 

ID If h can not be obtained, return Indication of 

fal lure. 

Ill) Compute, using the DELIMITER, c and d as In 
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c = LIMITCg(X), X, a, PLUS), 
d = LIMIT(gCX), X, b, MINUS). 



!v) Obtain f(Y) by assigning It a value computed 
from 

f(h(Y)) h»(Y). 

5. Solving Systems of Linear Algebraic Equations 

In sect. i|-5.2 the need to compute a polynomial Q(X) 
from a given one P(X) satisfying a given relation 

(1) QCX)-Q(X + 2 %P\ %\) » PCX) 

has been mentioned. 

The method of undetermined coefficients Is used to 

determine Q(X). Let Q(X) be a polynomial In X with degree 

one higher than P(X) and unknown coefficients CO, CI, ... , 

Cn. That Is 

n + 1 n 
Q(X) = X + Cn X + . . . + CI X + CO. 

Equation (1) will give n linear relations among these 
coefficients. Therefore the value of these C's can be 
obtained by solving the system of linear equations they must 
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satisfy. 

For the purpose of solving systems of linear algebraic 
equations/ A method known as the "Two-step f ractJon-f ree 
Gaussian elimination" [14] has been Implemented. This 
method Is an Improvement over a corresponding one-step 
method and features a procedure that keeps the size of the 
Intermediate expressions In the course of the reduction down 
by dividing them by a common factor which the procedure can 
predict. The advantage of this method over a more 
efficient and elaborate scheme [18] Is Its simplicity. 
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6. Convergence of Integrals 

There are many methods to determine the convergence or 
divergence of a given Improper tntegral. Some are limit 
tests/ others comparison tests etc. Some test for absolute 
convergence/ others conditional or uniform convergence. It 
Is a possible area for future work. There Is no real need 
to have such elaborate schemes in WANDERER, for after all it 
Is the value of the gtven Integral that Is desired. 
Cons Ider 




'a 

If f(X) Is a rational function of X, WANDERER combines 
convergence tests with evaluation algorithms as explained In 
chapters 4 and 5. If f Is not rational, then It Is not 
allowed to become Infinite except at a and b. This Is a 
convention on Inputs used by WANDERER. Thus, J Is a proper 
Integral If a, b, f(a) and f(b) are finite. If J Is 
improper, WANDERER has a test for absolute divergence. 

Test for absolute divergence: 

(1) If a and b are finite and b > a, then J Is 
absolutely divergent if 
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LlMIT((b-X)f(X),X,b, MINUS) t 0, 

or 

LIMlT((X-a)f(X),X,a,PLUS) t 0. 

(2) If a Is !nfln!te / then J Is absolutely divergent If 

LIMIT(X*f(X),X,a) t 0. 

(3) If b Is Infinite/ then J Is absolutely divergent If 

LlMlT(X*f(X) / X / b) f 0. 
If a function G(X) exists such that G'(X) = f<X), then 
the method of antlderl vatl ve can be used. That Is to com- 
pute J by evaluating 

LIMIT(G<X),X,b, MINUS) - LI Ml T(G(X) ,X,a, PLUS) . 

If this value Is finite, it Is the value of the given 
Improper Integral. If It does not exist, then J diverges. 
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CHAPTER VI I 
AN APPLICATION OF MACSYMA 
AND WANDERER 

0. Introduction 

The usefulness of a general purpose algebraic manipula- 
tion system such as MACSYMA In facilitating the solution of 
mathematical problems has been demonstrated. Using his 
"Symbolic Mathematical Laboratory" [15], Martin demonstrated 
solutions to three demanding problems In applied 
mathematics. These examples emphasize the fact that routine 
algebraic computation can be done by computer programs not 
only without error but much faster than by hand. The value 
of such a system Is especially appreciated when the expres- 
sions Involved are large and complicated. By employing 
such a computer facility the human problem solver may be 
freed from the tedtous and uninspiring manipulations to 
think more about the profound aspects of his problem. But 
this Is not all such a system can do. Moses's SIN [20]/ a 
program for Indefinite Integration/ provides a good example 
of successful mechanization of a mathematical process which 
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Is far from routine or straight-forward. The Introduction 
of SIN broadened the scope of algebraic manipulation systems 
slgnlf Icantl y. 

The purpose of this chapter Is to show how MACSYMA and 
WANDERER can be used to help solve complicated problems that 
are of practical Importance. One such problem Is the asymp- 
totic evaluation of certain contour Integrals arising In 
mathematical physics. Usually one starts with one or a set 
of differential equations describing a physical problem. 
Solving these equations by one method or another, most often 
by Integral transforms, one will arrive at a solution In the 
form of a definite Integral which Is often difficult. If not 
Impossible, to evaluate exactly. Frequently, one Is not so 
Interested In the exact solution but the behavior of the 
system when one parameter becomes very large or small. This 
is where asymptotic analysis is needed. 

MACSYMA Is used to obtain the asymptotic solution of an 
Infinite Integral. Many facilities provided In such an 
algebraic manipulation system can be Illustrated through 
this application. It Is also possible to show how the 
definite Integration capability provided by WANDERER Is 
needed for the successful computation of the results. 
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1. An Asymptotic Analysis Problem and an Outline of the 
Method of Steepest Descent 

Consider the infinite Integral 

-co 

C > 




where 



While Integrating J exactly may be Impossible/ Its 
asymptotic behavior as C becomes very large can be 
Investigated. To obtain the asymptotic expansion of J, the 
method of steepest descent [5,7] will be employed. 
Basically, the method of steepest descent consists in 
deforming the contour of Integration In such a way that the 
major contribution to the Integral arises from a small por- 
tion of the new path of Integration. The contribution will 
become more and more dominant as the parameter of Interest 
grows. This parameter here Is C. 

The first step In this method Is to find the new path 
of Integration. On a given contour, larger contributions 
comes from portions where the Integrand is larger in 
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absolute value and less oscillatory. Hence the requirements 
for a desirable path are: (a) the absolute value of the 
Integrand becomes maximum at a point, tl say, on the path, 
(b) The argument of the integrand is constant on the path 
near tl. The first requirement is obvious. The second Is 
essential, for if the, phase angle changes even slightly near 
tl, this change will be magnified by the very large factor C 
resulting in rapid oscillations of the Integrand and 
therefore negating any possible contribution from the point 
tl. Let U(tl,t2)= Rl (HC tl + % I t2)) and V(tl,t2) = lm(H(tl + %l 
t2)). Let p be a point where 



dU 
d tl 



dU 
d t2 



dV 
d tl 



dV 



d t2 



= 0. 



Then p is certainly a candidate for tl required In (a) and 
(b). Such a point is called a saddle point. The choice of 
the name saddle point will be made clearer later. The 
Cauchy-RIeman conditions Imply that criterion for such a 
point is M'(t)=0. There may be more than one such point in 
the complex t-plane. There are an infinite number of diffe- 
rent curves which pass through a saddle point and satisfy 
(a) and (b). Among them the path along which the integrand 
decrease in size most rapidly is the best. Tor the integral 
J, this means a curve on which U decreases most rapidly. 



Let H*(b) =0. It can be shown by use of properties of 
analytic functions that U varies most rapidly on curves 
lm(H(t))=const. (usually called level curves). If H"(b)ji0, 
two level curves will pass through the point b, Intersecting 
at right angle, as shown In fig. 1. 




fig. 1. 



On one of these two curves, say curve A, RUH(t)) Is minimum 
at b and Increases as t moves along A away from b. On the 
other hand, Rl(H(t)) Is maximum at t=b on the curve D and 
decreases as t moves away from b on D. Curve A Is called 
the steepest ascent path and D the steepest descent path. 
The point b Is usually referred to as a saddle point. If 
H"(b)j«0, b Is a saddle point of order 1, If H"(b)=0 and 
H IM (b)j«) of order 2, etc. If the original contour can be 
deformed onto one or a combination of such steepest descent 
paths, then the asymptotic expansion of the given Integral 
can be obtained by a rather routine procedure which 
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Involves: change of variable of Integration, Inversion of 
truncated power series at the saddle points and term by term 
Integrat Ion. 

Let us outline the steps of the solution procedure as 

fol lows. 

1) Locate and determine the order of saddle points 

of H(T). 

2) Compute U(tl,t2) and V(tl / t2) such that 
HCtl + %\ t2) = U(tl,t2) + SI V(tl,t2). 

3) Obtain V(tl,t2) - const, curves .which pass 
through the relevant saddle points. 

k) Examine the V(tl/t2) « const, curves to 
determine whether deformation of contour can be made 
to curves through the saddle points. 

5) Change the variable of Integration. 

6) Express t as a truncated series In the new 
variable about each relevant saddle point. 

7) Determine the coefficients In the above series. 

8) Apply Watson's lemma to obtain the first few 
terms of the asymptotic expansion by Integrating 
term by term. 
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2. Sol utlon Steps 

Presented here are the solution steps of the asymptotic 
expansion of iKC) In the exact sequence as they have been 
carried out using MACSYMA. The lines labelled (CI) are 
Input commands. A command line ends wtth either a @ sign or 
a $ sign. The @ sign causes results obtained by executing 
the command line to be displayed in a subsequent line 
labelled (Dl). A (CD line together with a corresponds 6DJ) 
line will be referred to as step I. The $ as an end of 
command line character suppresses display of results for 
that line. Explanatory texts will be Inserted between 
lines. To avoid becoming a user's manual for MACSYMA, 
explanation for the commands used In the solution will be 
made quite brief. For a more detailed look at MACSYMA the 
reader Is referred to [16], 

(CI) P:8$ 

P Is a parameter which Is set depending on the number 
of terms desired in the asymptotic expansion. By setting P 
to 8, we shall obtain the first 4 terms. The reason will 
become evident later. 

(C2) H(T):=-T**4/4-2l*T@ 
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T 

(D2) H(T): %\ T 

k 

defining the function H(T) 

(C3) DIFF(H(T),T)=0@ 

3 
(D3) - T - %\ = 

creating an equation H'(T)=0 



(CU) SOLVE(%,T)Q 
SOLUTION 



SQRTC3) - %\ 
(E4) T = 



2 

(E5) T = %\ 

- SQRT(3) - %\ 
(E6) T = ■ 



2 

(D6) (EU,E5,E6) 

The £ sign used In (CU stands for the last (03), In 
general a % sign represents the last expression labelled 
(Dl). SOLVE In (CU) Is an Invocation of the MACSYMA SOLVE 
program (see [10] for Its capabilities and limitations). 
The roots of H'(T)=0 give three first order saddle points. 
Passing through each of these points there will be one 
steepest ascent and one steepest descent level curve 
(Im (H(t))= const, curves). The saddle points In (E'O and 
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(E6) are the points B and A respectively shown In fig. 2. 
The value of H at these three saddle points will now be 
computed. 

(C7) HB:RATSIMP(H(PART(EI*,2)>>@ 

- 3 %\ SQRTC3) - 3 

(D7) 

8 

What has been done In step 7 Is the computation and 
simplification of H(B). The command RATSIMP causes rational 
simplification which Is essentially putting expressions to 
be simplified Into the form of one numerator. and one 
denominator and perform all possible GCD cancellations. The 
command PART allows a user to obtain subexpressions of an 
expression. PART(E^/2) returns the second part of the 
equation EU which Is I ts' right-hand side. Note that 
commands can be nested. 

(C8) HC:RATS!MP<H(PART(E5,2)))@ 

3 
(D8) 

(C9) HA:RATSIMP(H(PART(E6,2)))@ 

3 %\ SQRTC3) - 3 



(D9) 

8 

The point T=SI turns out to be Irrelevant because the 
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new path which will be determined does not pass through it. 
The next goal Is to obtain curves passing through the saddle 
points (SQRT(3)/2 / -1/2) and (-SQRT(3)/2 / -1/2) along which 
the Imaginary part of H(T) Is constant and the real part of 
H(T) varies most rapidly/ I.e. the steepest paths through 
the saddle points. 

(CIO) EXP:EXPAND(H(T1+$I*T2))@ 

k 2 2 

T2 3 3 Tl T2 

(D10) ♦ %\ Tl T2 + - - 

k 2 

3 Tl 
-%\ Tl T2 + T2 %\ Tl 

k 

(Cll) V(T1/T2):="COEFF(D10/%I)@ 

3 3 
(DID V(T1 / T2):= Tl T2 - Tl T2 - Tl 

By steps 10 and 11/ the Imaginary part of H(T1+ %\ T2) 
Is found and given a function name V(T1/T2). The command 
COEFF(exp/Var/n) computes the coefficient of var**n In exp. 
In thts case, the coefficient of %\ In (D10) Is exactly the 
Imaginary part of (D10). 

(C12) U(T1/T2):="COEFF(EXP/2;I/0)@ 
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4 2 2 4 

T2 3 Tl T2 Tl 

(D12) U(T1,T2): = + + T2 

4 2 4 

U(T1,T2) ts defined to be the real part of 
H(T1 + %\ T2). 



(C13)V(Tl,T2)=V(SQRT(3)/2, -1/2 )@ 

3 3 3 SQRT(3) 
(D13) Tl T2 - Tl T2 - Tl = 

8 



Obtained In (D13) Is the equation of a curve, V = 
constant, which passes through the saddle point B, 
(SQRT(3)/2 / -1/2). This curve will be referred to as CR. 

(CHOV(Tl,T2)=V(-SQRT(3>/2, -1/2 )@ 

3 3 3 SQRT(3) 
(D14) Tl T2 - Tl T2 - Tl = 

8 

(D14) Is the equation of a curve, V=constant, pass tng 
through the saddle point A, (-SQRT(3)/2 / -1/2). Let this 
curve be CL. The curves CL and CR have to be examined 
carefully by the human problem solver to determlng the new 
path of Integration. The manner In which they extend to 
Infinity ts often Important In deforming the contour. 
Asymptotes to these curves can be found easily. For an 
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It can be seen In ftg. 2 that by use of Cauchy's 
Integral theorem/ the original contour can be deformed to 
the contour Ca+Cb. For the purpose of change of path of 
integration a new variable R ts Introduced which will be the 
parameter of our steepest descent paths. 

(C17) -R**2=H(T)-H(T0)Q 

k 
2 T 1 

(D17) " R - *l TO %\ T + - 

TO stands for either of the two saddle points A or B. 
Indeed if T0=A, R ts real If and only If T Is a point on Ca, 
for lm(H(T)-H(A))=0 and ss (H(T)-H(A)) only for points on 
Cj.. The same can be said about T0=B. T has to be expressed 
In terms of R In order to perform the change of variable. 
Solving for T as a function of R exactly In (D17) Is not 
necessary. What Is needed Is the first few terms of a power 
series expansion of T In terms of R about the point T=T0 
(the saddle point). 

(C18)T=T0+DOSUM(J,l,P/B[J]*R**J)@ 

8 7 6 5 h 

(D18)T =T0+BR+BR+BR+BR+BR 

8 7 6 5 k 

3 2 

+ B R +BR +BR 
3 2 1 
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T Is set to a truncated power series In R. The unknown 

coefficients Bl will be solved for by use of (D17). 

(C19) E:SUBSTITUTE(D18 / R**2+PART(D17, 2) )@ 

8 7 6 5 1* 
(D19) -(TO+B R +B R +B R + B R +B R 

8 7 6 5 1* 

3 2 1+ 8 7 

+ B R + B R + B R) A - %l (TO + B R + B R 
3 2 1 8 7 

6 5 4 3 2 

+ BR+BR+BR+BR+BR+BR) 
6 5 4 3 2 1 

2 1 
+ %\ TO ♦ R + - 

k 

The expression (D19) has been given a name E. Before E 

Is expanded In order to collect terms, a few simplification 

rules are defined so as to discard powers of R higher than 8 

In the expansion process. This greatly reduces the 

Intermediate expression bulge which would otherwise occur. 

(C20) DECLARE (N / MPRED)$ 

(C21) MPRED(X):=IF X > P THEN TRUE ELSE FALSE$ 

(C22) TELLSIMP (R**N,0)$ 



IS THE REPLACEMENT 

A simplification rule has been set up so that any R**n 
Is replaced by If N > p which Is 8. 
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(C23) TELLSIMP(T0**3,-%I )$ 

An additional simplification rule Is defined so that 
T0**3 sha'il be replaced by -%\ automatically. 

(C2U) RATVARS(R)$ 

R shall be the main variable In subsequent rational 
s impl i f I ca 1 1 ons . 



(C25) E:PART(RAT(E),1)Q 
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(D25) 

8 2 2 

R((-12B B -12B B -12B B -6B )T0 
17 2 6 3 5 k 

2 2 

+ ( - 12 B B - 24 B B B +C-2UB B -12B )B 
16 125 13 2 4 

2 3 2 2 2 

■ 12 B B ) TO - U B-12B BB-6B B 
23 15 12 4 13 

2 4 

- 12 B B B - B ) 
12 3 2 

7 2 

+ R (( - 12 B B - 12 B B - 12 B B ) TO 
16 2 5 3 4 

2 2 2 

+ ( - 12 B B - 24 B B B - 12 B B - 12 B B ) TO 

1 5 1 2 4 13 2 3 

3 2 3 

■H B - 12 B B B - 4 B B ) 
1 4 12 3 12 

6 2 2 

+ R (( - 12 B B - 12 B B ■ 6 B ) TO 
15 2 4 3 

2 3 

+ ( - 12 B B ■ 2li B B B - li B ) TO 
1 4 12 3 2 

3 2 2 

- 1> B B - 6 B B ) 
13 12 



5 2 

+ R (( - 12 B B - 12 B B ) TO 
14 2 3 
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2 2 3 

+ ( - 12 B B-12BB)T0-UB B) 
13 12 12 



k 2 2 2 4 

+R (( - 12 B B -6B )T0 - 12 B B TO-B ) 
13 2 12 1 

3 2 3 2 2 2 

+ R ( - 12 B B TO - 4 B TO) + R (It - 6 B TO ) 

1*2 1 1 



E Is set to the huge expression above which Is the 
numerator of a truncated expansion of E; the denominator Is 
a constant. 

(C26)FOR J:l STEP 1 UNTIL J > (P-l) 

DO BLOCK (EQ[J]:COEFF(E,R,J+1)=0,DISPLAY(EQ[J]>>$ 

This do loop Is used to generated the coefficients of 
the various powers of R In E and store them In the form of 
equations In an array RL. Each EQ[J] will be displayed 
after It Is set. 



2 2 

k - 6 B TO = 
1 

2 3 

- 12 B B TO - k B TO - 
12 1 

2 2 2 k 

( - 12 B B - 6 3 ) TO - 12 B B TO - B = 
13 2 12 1 
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( - 12 B B - 12 B B ) TO 
1^23 

2 2 3 

♦ ( - 12 B B - 12 B B ) TO - It B B 
13 12 12 

2 2 

(-12B B -12B B -6B )TO 
15 2 k 3 

2 3 

+ ( - 12 B B - 21 B B B - k B ) TO 
11* 12 3 2 

3 2 2 

- k B B-6B B » 
13 12 



( - 12 B B - 12 B B - 12 B B ) TO 
16 2 5 3 k 

2 2 2 

+(-12B B-24BBB-12BB - 12 B B) 
15 121 13 23 

3 2 3 

TO ■ l( B B-12B BB-UBB =0 
1 it 12 3 12 



2 2 

( - 12 B B -12B B -12B B -6B )T0 
17 2 6 3 5 k 

2 2 

+ ( - 12 B B - 24 B B B +C-21B B -12B ) 
16 125 13 2 

2 3 2 

B - 12 B B ) TO - U B-12B BB 
*» 23 15 12 '♦ 
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2 2 2 k 

- 6 B B -12 B B B-B =0 
13 12 3 2 



Being conscious of the lack of storage space left for 

subsequent manipulations, we get rid of the computation 

history to create some space. 

(C27) KILL(HISTORY)@ 

(CI) S0LVE(EQ[1],B[1])@ 

SQRT(6) 

(El) B = 

1 3 TO 

SQRT(6)« 



(E2) B - 

1 3 TO 

(D2) (E1,E2) 

The choice of the value for Bl here affects only the 

sense with which the steepest decent paths are traversed. 

Thus If we choose (El), R would vary from IMF to MINF on Ca 

and from MINF to INF on Cb. 

(C3) B[1]:PART(E1,2)$ 

A do loop Is used In (CO to solve for the remaining 
B's (B2 through B7). 
(CO FOR J:2 STEP 1 UNTIL J > (P-l) 

DO BLOCK (SOL:EV(SOLVE(EQ[J] / B[J]) / EVAL) / 
B[J] :PART(S0L,2),DISPLAY(S0L))$ 

The EV command, with EVAL as the third argument, causes 
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The unknown coefficients B! have been determined. The 
next goal Is to compute 
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r™ 2 dT 

' EXP(-C R ) -- dR 



dR 
-co 

dT 
for TO^A and TQ = B„ Odd terms of R In -- do not contribute. 

dR 

Therefore only k terms need be Integrated. These terms are 

even In R, thus the range of Integration can be changed to 

(0 IMF). A function FN(X) will be defined for carrying out 

this term-by- term integration at a variable point T0=X. 

(C5) FN(X) :=BL0CK(AN5:C, FOR J:l STEP 2 UNTIL J > (P-l) 

DO ANS:AN5+J*RATSIMP(EV(B | J |,T0=X))* 

DEFINTUE**(-C*R**2)*R**(J-1),R,0, I NF ), RETURN (ANS )) $ 

(CG) FN(~SQRT(3)/2-^l/2)@ 

IS THE EXPRESSION 

- C 
POSITIVE, NEGATIVE, OR ZERO 

NEGATIVES 

For the purpose of integration, WANDERER askes the user 

about the sign of C. The contribution from saddle point A 

Is given in the next result. 



SQRT(G) SQRTUPI) 
( DG) _ ____„ 

(3 SQRT(3) + 3 ?jl ) SQRT(C) 
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3 SQRTUPI) (3 SQRTC3) + 23 %\) 

+ 3/2 

J» SQRT(6> (27 %\ SQRT(3) + 27) C 

5 SQRTUPI) (25 SQRT(3) + 247 %\) 

5/2 
5184 SQRT(6) C 

- (35 SQRTUPI) (SQRT(6) (672 SQRT(3) + 1816 %\) 

3582 SQRT(3) 13077 %\ 7/2 

))/(839808 C ) 

SQRT(6) SQRT(6) 



(C7) FN(SQRT(3)/2-%!/2)@ 

IS THE EXPRESSION 

- C 
POSITIVE, NEGATIVE, OR ZERO 

NEGATIVEQ 

SQRT(6) SQRTUPI) 
(D7) - 

(3 SQRT(3) -331) SQRT(C) 



3 SQRTUPI) (3 SQRT(3) - 23 %\) 

* 3/2 

4 SQRT(6) (27 %\ SQRT(3) - 27) C 

5 SQRTUPI) (25 SQRT(3) - 247 %\) 

+ _«___-_ — »_«_____ 

5/2 
5184 SQRT(6) C 



+ (35 SQRTUPI) (SQRT(6) (672 SQRT(3) - 1816 %\) 

3582 SQRT(3) 13077 %\ 7/2 

+ --))/(839808 C ) 

SQRT(6) SQRT(6) 
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The above Is the contribution from saddle point B. Now 

the final result Is obtained by computing the following 

(-3 /3 % I -3>/8 (3 /3 %l-3)/8 

2 UE D7 - 2 %E D6) 

where D6 and D7 are as given above. This expression Is, by 
Inspection, equivalent to 

(-3 >/3 % I - 3 > / 8 (3 /3 *l-3)/8 _ 
2 UE D7 + %E D7) 

where D7 Is the complex conjugate of D7. Therefore/ 

(C8) i»*RL(%E**HA*D7)Q 

3 SQRT(3) %\ - 3 

8 SQRT(6) SQRT(%PI-> 
(D8) k RLUE ( - - 

(3 SQRT(3) - 3 %\) SQRT(C) 

3 SQRTUPI) (3 SQRT(3) - 23 %\) 

+ _«__--__>-.«. «.„ ____-_-_--- 

3/2 

k SQRT(6) (27 %\ SQRT(3) - 27) C. 

5 SQRTUPI) (25 SQRT(3) - 247 %\) 

5/2 
5184 SQRT(6) C 

+ (35 SQRTUPI) (SQRT(6) (672 SQRT(3) - 1816 %\) 

3582 SQRT(3) 13077 %\ 111 

+ — ))/(839808 C )) 

SQRT(6) SQRT(6) 
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This Is the first four terms of the desired asymptotic 
expansion. RL is not a command of MACSYMA. It is used here 
to denote the real part of an expression. 
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CHAPTER VI I I 
CONCLUSIONS AND SUGGESTIONS FOR FURTHER WORK 

The evaluation of definite Integrals Is a classical 
problem In mathematics. Great Ingenuity Is frequently 
requlred / with many Integrals demanding special devices. 
The lack of a sufficiently general theory makes evaluation 
of definite Integrals very difficult. It Is doubtful that a 
theory/ comparable In generality to the Rlsch Integration 
algorithm [26] for Indefinite Integrals, can be developed In 
the near future. We have shown that the convergence of a 
class of Integrals of elementary functions Is recursively 
undecldable. The proof Is In appendix E. The WANDERER pro- 
gram presented here Is a prototype heuristic computer pro- 
gram for the symbolic evaluation of definite Integrals. It 
contains both general methods such as contour Integration, 
residue theory and differentiation with respect to a 
parameter, In addition to quite a few special methods for 
specific types of Integrals. Clues as to which method to 
use for a given problem are obtained from the Integration 
range and the form of the Integrand. Although many types of 
definite Integrals can be evaluated by WANDERER, It, as 
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almost any other heuristic computer program, has Its limita- 
tions. Foremost among these Is the fact that It Is a pro- 
gram for the evaluation of real definite Integrals of 
elementary functions. 

It Is hoped that the work reported here may provide a 
starting point for new approaches to the evaluation of 
definite Integrals from the viewpoint of symbol manipula- 
tion. The advantage of this approach Is twofold: (1) A com- 
puter can use integration methods that are too lengthy or 
complicated to be carried out by hand; (2) Such a computer 
program contains a collection of powerful methods that can 
interact with one another and can produce answers to Inte- 
grals not present in any finite table. For this reason, a 
good definite integration program together with other 
facilities provided In an algebraic manipulation system can 
be very useful to applied mathematicians. 

DELIMITER is a rather sophisticated program for comput- 
ing limits of elementary functions. Such a program has been 
shown to be a very useful tool in an algebraic manipulation 
system. The method of comparing orders of Infinity used In 
DELIMITER Is an Important concept which Is useful In places 
other than the computation of limits. 
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A basic assumption of DELIMITER Is that the operator 
LIMIT Is distributive over the operators + , -, *, and **, 
I.e., rules such as the limit of a sum Is the sum of the 
limits hold. These rules are valid only if the subproblems 
thus generated produce answers which do not lead to an 
Indeterminate form. Of course there Is L'Hospital's rule 
which can be applied to some Indeterminate forms, but 
generally only those Involving only the operator * or **. 
Although many aspects of the Indeterminate form (INF-INF) 
have been considered. In chapter 2, the algorithm for Its 
determination Is not complete. A powerful series expansion 
program would be helpful In some cases. Yet, such a pro- 
gram can not solve all problems. Consider for Instance, the 
1 Im! t problem 

2 2 
LIMIKSIN X + COS X,X,INF). 

To obtain such limits, an algebraic manipulation system 
must be able to detect all constant Identities. This Is not 
possible for the set of all expressions [22]. It may be 
possible for a proper subset of all expressions. For 
example, many trigonometrical Identities disappear If all 
trigonometric functions are radically transformed Into sums 
of complex exponentials. This Is a basic problem of great 
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practical and theoretical Interest [61. Such simplification 
capabilities would be of use to many other programs as well. 

DELIMITER can be augmented by Increasing the types of 
functions It can handle which may Include functions defined 
by integrals. Another possible area of research is the 
automatic determination of superior and inferior limits. 

The fact that WANDERER can evaluate many non-trivial 
definite Integrals does not mean that It can compete with an 
expert human Integrator yet. For one thing, a mathematician 
can usually construct a function of a complex variable and a 
suitable contour for evaluating different Integrals. 
WANDERER cannot form a contour based on analysis of a given 
real Integral. It simply selects from the cases known to 
It. WANDERER would be much more powerful If It could 
determine, for a given Integral, whether the method of con- 
tour integration and residue theory were applicable, and, If 
It were, evaluate the Integral by forming a suitable Inte- 
gral around a closed contour. 

The evaluation of contour Integrals by residue theory 
usually requires the solution of algebraic or transcendental 
equations. WANDERER uses the SOLVE program of MACSYMA for 
this purpose. SOLVE has Its limitations and Its Improvement 
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is also a possible area for further work. 

Suppose one wants to evaluate a closed contour Integral 
by residue theory. An Interesting general question is: what 
knowledge about the integrand or the functions used in form- 
ing it Is necessary. We think the following are essential: 

( 1 ) eva 1 uat ion; 

(2) differentiability and derivatives; 
( 3 ) si ngul ar i t i es; 

CO asymptotic behaviour. 

A natural extension of the work here is integration 
over arbitrary user-specified contours. This should not be 
difficult to do, except for the lack of notation. 
Specifically, we must specify an arbitrary contour to a com- 
puter and devise a general data structure for use In 
representing contours. For Instance, consider the 
specification of the following indented contour. 
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One obvious way is to specify a contour piecewlse. 
Each piece of curve has a parametric form and a starting and 
ending value for the parameter. There may be other 
approaches. 

Further work In this area of symbol manipulation might 
tncl ude: 

(1) Design of computer algorithms for testing convergence 
and divergence of Integrals/ 

(2) Summation of Infinite series by residue computations/ 

(3) Investigation of algorithms for deflntte. Integration of 
special functions. 

(h) Evaluation of multtple definite Integrals. 

(5) Application of definite Integration programs In solution 

of differential equations. 
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APPENDIX A 
A FLOWCHART OF DELIMITER 

A flowchart Is presented In the next page which details 
the flow of control of DELIMITER. Listings of the programs/ 
written In the LIS? programming language/ may be obtained 
from the author. 

The routine LIM/ appearing In the flowchart/ Is a 
program which applies the 'trivial' limit rules/ makes use 
of subroutines to compute limits of sums/ products/ powers, 
and the functions SIN/ COS/ TAN/ LOG/ SINH/ COSH/ TANH. It 
calls LIMIT recursively and makes use of L'Hospltal's rule 
and other routines when needed. 
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LIMIT(E(X),X,L) 
rSTART^ 




L'HOSPITAL'S RULE 
FAILED? 



HRETURN ANSWERI 



L'HOSPITAL'S RULE 



COMPARE ORDERS 
OF INFINITY OF 
N AND D 

1 



N AND D CONTAIN 
FUNCTIONS IN THE 
FORM e**R(X) 
WHICH TEND TO INF? 



E RATIONAL 
IN X ? 



RATIONAL 
FUNCTION 
ROUTINE 



N-NUMERATOR 
OF E 
D-DENOMINATOR 
Of E 



RETURN LIM(NL/DL) 



Tf 



NLAMPDL= INFINITY? \ 



I D = 1 ?F 



RETURN 
LIM(N) 



N AND D ARE 
RP-FQRMS? 



RP-FORM ROUTINE 



V 



NL=INF OR MINF 

AND ? 
DL=INF OR MINF 



DL^S} 3 ""*® 



NL=0 AMD DL«0? 



NL-LIMIT(N,X,L) 
DL-LIM1T(D,X,L) 






D— + ?l 



D-*0- 71 



HL>0 



1NL>0 ?f 



3- 



RETURN IN3 



RETURN MINFT - 



LRETURN INFINITY 
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APPENDIX B 
A PROOF 



Infinite Integrals of 



n 



%\ k Z m 
f(Z) « %E Z 

where n > 0, k f real, Rl (m) > -1, and n - RICm) > 1, 

have been dlscused In chapter U-3.2. Some results derived 

there depends on the proof of (1) and (2) below. 

The objective here Is to supply the proof of 

1) limit / f(Z) dZ = 0, 



and 



limit f 
fmlt f 



2) limit / f(Z) dZ = 0. 



where CR and C are the circular contours given below, 



->x 
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Let 



Z = r e and m = a + lb. 



then 



\L 



f(Z) dZ 



(i) 



m+l / n Iile 

! r / EXPCt kr e + i f nt + i) * ) de 

a +l f"X/Z7l n 

^ r / EXP( - k SIN(n e) r - b e) de 



= M 



For ss e «S %P I /( 2 n) we have 



n 



Thus 



EXP( - k SINCn e) r ) ^ 1. 



rit/zn 
a+l/ -be 

M s r / e de 



Since a + 1 > 0/ this completes the proof of (1). 



Now from equation (I) we have 
a+l 



•tt/2. 



n be 



M = 



n 



EXP( - k SIN(e) r ) de 

n 



a+l 



-7C/2 



2 k e r be 



EXP( - 



) de 



%P\ 
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a + 1 



5rP I r 



n SPI b 



■11 - EXP( - k r 



■)) 



2 n 



The fact that 



k n r + % P i b 

1 > completes the proof of (2), 
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APPENDIX C 

PERFORMANCE OF WANDERER 
AND DELIMITER 

As a measure of the performance of WANDERER and the 
limit programs, selected problems have been timed. The time 
sharing system used is the ITS of the Artificial 
Intelligence Laboratory at MIT which uses a PDP-10 computer 
with a memory cycle time of about 2.75 microseconds. Time 
used for parsing the input string and display of the 
computed result has been excluded In order to obtain an 
approximation to the time actually spent Inside WANDERER or 
DELIMITER. Garbage collection In the LISP [17] system In 
which MACSYMA is written may take place during a 
computation. Although garbage collection Is a slow 
procedure/ It Is only fair to regard It as part of the 
computation process being carried out. Therefore/ the time 
it required has been included In the timing experiments. 
The results are put In the form of two tables. An * Is used 
in the tables to indicate computation requiring LISP garbage 
col 1 ect Ion. 
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TABLE I 

-PRFORMANCE OF DELIMITER 



L. i n i t 



puted 



Rosu i t 
e q u I v a I o fit Co 



Tine I r 1 
Sec , 



LIMIT X 
X-*Q + 



0,32 



LIMIT ( 1 + .i / x ; 



2E 



1.U2 



IMI T 

-»+ CO 



)0(X) ♦ %£ 

sqp.TCx + i) 



- \ 



) + °<E 



3,00 



LIMIT IE 
X-*+ao 



+ I J X 
- %E 



NF 13.18 



LIMIT 



- G 



NF 



1.15 



LIMIT - - - 
X-ff/ 2 C 



1. 52 






■ 2L 

2 



;qr- 



N ! F 



Integral computed 



fr 



-CD 



dX 



2 
'° X + 3 X * 2 



/ 



"* dX *PI 



20 20 SJNUPI/20) 

° X +1 
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TABLE 1 1 




PERFORMANCE OF WANDERER 




Resul t 


Ttme fn 


equivalent to 


Sec. 



dX (principal) 0.55 



L0G(2) 2.02 



, 2 

X + A X + B %?\ B + 3 %P I * 

• dX 7.1*2 

hi 12 

'-"X + 10 X +9 



2 
*X + A X + B 3 A LOG(3) + (B +.3)%PI 

dX - 11.32 

l* 2 2U 



* 



X + 10 X +9 



2 

COS (X) - SIN(X) dX %PI 1.67 



1.08 



X +1*9 
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-co 



3 D -5 SQRT(X) 
X °^E dX 



D 
(2016 |E ) / 78125 



0.95 



.oo 2 

COS (X) 



dX 



% P I / 2 8 



1.39 



-oo 



dX 



SQRT(X) (2 X + 1) 



%?\ /SQRT(2) 



0.81 



-CD 



1/3 

X dX 



(It SQRT(X) + 3) 



3/2 



DIVERGENT 



6.62 
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APPENDIX E 
AN UNDECIDABILITY RESULT 

Let SI = {P(X1,X2, .. .,Xn)} be a set of polynomials with 
Integral coefficients In XI, ... , Xn and S2 a set of 
functions F of the form 

F(X1, ... , Xn) 

2 2 £ 2 2 

= (n+1) (P (XI, ... , Xn) + £ (SINrtXj)KJ (XI, ... , Xn))-1, 

H 

where Kj Is the dominating function |6| for 



d 2 

P. 

dXj 



Richardson has shown |22,6| that 

Lemma 1. For F In S2 the predicate "there exists an n-tuple 
B of real numbers such that F(B) < 0" is recursively 
undec Idabl e. 

Lemma 2. If F(B) < for some n-tuple B real numbers then 
there exists an n-tuple A of nonnegatlve Integers such that 
P(A)=0 and therefore F(A) = -1. 
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Lemma 3, Let 



h(X) - X SIN(X) and g(X) = X SiN(X ) 

Then for any real numbers Al, ... , An and any < 6 < 1 

there exists b > such that 

|h(b)-Al|<€, |h( E (b))-A2|<6 / ... , | hC g( ... g( b )...)) -An |< €. 

By use of these Lemmas we can show 

Corollary. For any G in a set of functions of the form 

G(X) - F(h(X),h(g(X)), ... , h(g(.. .(g(X))...))) + 1/2, 
the predicate "there exists a real number t such that G(t) = 
0" is recursively undecidable. 

Proof : Suppose the predicate is recursively decidable, then 
we have 

1) If there exists a t such that G(t)=0 then there exist 
real numbers Al, ... , An such that F(A1, ... , An) < 0; 

2) If there exist real numbers Al, ... , An such that 
F(A1, ... , An) < then (by lemma 2) there exist 

n on negative integers Bl, ... , Bn such that 

F(B1, ... , Bn)=-1. Thus, there exists a real number c such 
that G(c) < 0. Since G(X) can be large and positive and 
G(X) is continuous, this Implies that there exist a t such 
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that G(t)=0. 

This contradicts lemma 1. 

From the definition of the Function F, one can see G is 
always -1/2 and G Is large and positive except at the 
vicinity of a finite number of points where G Is negative. 

Theorem ! The convergence of a set of integrals of the form 

■°> dX 



2 2 
'-a>(X +1)G (X) 



Is recursively undecidable. 

Proof: This Integral Is convergent If and only If G(X) has 
no real zero. 
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